Index: sansmodels/src/sans/models/CoreShellCylinderModel.py
===================================================================
--- sansmodels/src/sans/models/CoreShellCylinderModel.py	(revision ae3ce4eaa593202ebf49122e196897c0d45a25b6)
+++ sansmodels/src/sans/models/CoreShellCylinderModel.py	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
@@ -1,11 +1,23 @@
 #!/usr/bin/env python
+"""
+	This software was developed by the University of Tennessee as part of the
+	Distributed Data Analysis of Neutron Scattering Experiments (DANSE)
+	project funded by the US National Science Foundation.
+
+	If you use DANSE applications to do scientific research that leads to
+	publication, we ask that you acknowledge the use of the software with the
+	following sentence:
+
+	"This work benefited from DANSE software developed under NSF award DMR-0520547."
+
+	copyright 2008, University of Tennessee
+"""
+
 """ Provide functionality for a C extension model
 
 	WARNING: THIS FILE WAS GENERATED BY WRAPPERGENERATOR.PY
- 	         DO NOT MODIFY THIS FILE, MODIFY core_shell_cylinder.h
+ 	         DO NOT MODIFY THIS FILE, MODIFY ../c_extensions/core_shell_cylinder.h
  	         AND RE-RUN THE GENERATOR SCRIPT
 
-    @author: Mathieu Doucet / UTK
-    @contact: mathieu.doucet@nist.gov
 """
 
@@ -16,5 +28,5 @@
 class CoreShellCylinderModel(CCoreShellCylinderModel, BaseComponent):
     """ Class that evaluates a CoreShellCylinderModel model. 
-    	This file was auto-generated from core_shell_cylinder.h.
+    	This file was auto-generated from ../c_extensions/core_shell_cylinder.h.
     	Refer to that file and the structure it contains
     	for details of the model.
@@ -78,4 +90,13 @@
         
         return CCoreShellCylinderModel.runXY(self, x)
+        
+    def set_dispersion(self, parameter, dispersion):
+        """
+            Set the dispersion object for a model parameter
+            @param parameter: name of the parameter [string]
+            @dispersion: dispersion object of type DispersionModel
+        """
+        return CCoreShellCylinderModel.set_dispersion(self, parameter, dispersion.cdisp)
+        
    
 # End of file
Index: sansmodels/src/sans/models/CoreShellModel.py
===================================================================
--- sansmodels/src/sans/models/CoreShellModel.py	(revision ae3ce4eaa593202ebf49122e196897c0d45a25b6)
+++ sansmodels/src/sans/models/CoreShellModel.py	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
@@ -1,11 +1,23 @@
 #!/usr/bin/env python
+"""
+	This software was developed by the University of Tennessee as part of the
+	Distributed Data Analysis of Neutron Scattering Experiments (DANSE)
+	project funded by the US National Science Foundation.
+
+	If you use DANSE applications to do scientific research that leads to
+	publication, we ask that you acknowledge the use of the software with the
+	following sentence:
+
+	"This work benefited from DANSE software developed under NSF award DMR-0520547."
+
+	copyright 2008, University of Tennessee
+"""
+
 """ Provide functionality for a C extension model
 
 	WARNING: THIS FILE WAS GENERATED BY WRAPPERGENERATOR.PY
- 	         DO NOT MODIFY THIS FILE, MODIFY core_shell.h
+ 	         DO NOT MODIFY THIS FILE, MODIFY ../c_extensions/core_shell.h
  	         AND RE-RUN THE GENERATOR SCRIPT
 
-    @author: Mathieu Doucet / UTK
-    @contact: mathieu.doucet@nist.gov
 """
 
@@ -16,5 +28,5 @@
 class CoreShellModel(CCoreShellModel, BaseComponent):
     """ Class that evaluates a CoreShellModel model. 
-    	This file was auto-generated from core_shell.h.
+    	This file was auto-generated from ../c_extensions/core_shell.h.
     	Refer to that file and the structure it contains
     	for details of the model.
@@ -72,4 +84,13 @@
         
         return CCoreShellModel.runXY(self, x)
+        
+    def set_dispersion(self, parameter, dispersion):
+        """
+            Set the dispersion object for a model parameter
+            @param parameter: name of the parameter [string]
+            @dispersion: dispersion object of type DispersionModel
+        """
+        return CCoreShellModel.set_dispersion(self, parameter, dispersion.cdisp)
+        
    
 # End of file
Index: sansmodels/src/sans/models/EllipsoidModel.py
===================================================================
--- sansmodels/src/sans/models/EllipsoidModel.py	(revision ae3ce4eaa593202ebf49122e196897c0d45a25b6)
+++ sansmodels/src/sans/models/EllipsoidModel.py	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
@@ -1,11 +1,23 @@
 #!/usr/bin/env python
+"""
+	This software was developed by the University of Tennessee as part of the
+	Distributed Data Analysis of Neutron Scattering Experiments (DANSE)
+	project funded by the US National Science Foundation.
+
+	If you use DANSE applications to do scientific research that leads to
+	publication, we ask that you acknowledge the use of the software with the
+	following sentence:
+
+	"This work benefited from DANSE software developed under NSF award DMR-0520547."
+
+	copyright 2008, University of Tennessee
+"""
+
 """ Provide functionality for a C extension model
 
 	WARNING: THIS FILE WAS GENERATED BY WRAPPERGENERATOR.PY
- 	         DO NOT MODIFY THIS FILE, MODIFY ellipsoid.h
+ 	         DO NOT MODIFY THIS FILE, MODIFY ../c_extensions/ellipsoid.h
  	         AND RE-RUN THE GENERATOR SCRIPT
 
-    @author: Mathieu Doucet / UTK
-    @contact: mathieu.doucet@nist.gov
 """
 
@@ -16,5 +28,5 @@
 class EllipsoidModel(CEllipsoidModel, BaseComponent):
     """ Class that evaluates a EllipsoidModel model. 
-    	This file was auto-generated from ellipsoid.h.
+    	This file was auto-generated from ../c_extensions/ellipsoid.h.
     	Refer to that file and the structure it contains
     	for details of the model.
@@ -72,4 +84,13 @@
         
         return CEllipsoidModel.runXY(self, x)
+        
+    def set_dispersion(self, parameter, dispersion):
+        """
+            Set the dispersion object for a model parameter
+            @param parameter: name of the parameter [string]
+            @dispersion: dispersion object of type DispersionModel
+        """
+        return CEllipsoidModel.set_dispersion(self, parameter, dispersion.cdisp)
+        
    
 # End of file
Index: sansmodels/src/sans/models/EllipticalCylinderModel.py
===================================================================
--- sansmodels/src/sans/models/EllipticalCylinderModel.py	(revision ae3ce4eaa593202ebf49122e196897c0d45a25b6)
+++ sansmodels/src/sans/models/EllipticalCylinderModel.py	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
@@ -1,11 +1,23 @@
 #!/usr/bin/env python
+"""
+	This software was developed by the University of Tennessee as part of the
+	Distributed Data Analysis of Neutron Scattering Experiments (DANSE)
+	project funded by the US National Science Foundation.
+
+	If you use DANSE applications to do scientific research that leads to
+	publication, we ask that you acknowledge the use of the software with the
+	following sentence:
+
+	"This work benefited from DANSE software developed under NSF award DMR-0520547."
+
+	copyright 2008, University of Tennessee
+"""
+
 """ Provide functionality for a C extension model
 
 	WARNING: THIS FILE WAS GENERATED BY WRAPPERGENERATOR.PY
- 	         DO NOT MODIFY THIS FILE, MODIFY elliptical_cylinder.h
+ 	         DO NOT MODIFY THIS FILE, MODIFY ../c_extensions/elliptical_cylinder.h
  	         AND RE-RUN THE GENERATOR SCRIPT
 
-    @author: Mathieu Doucet / UTK
-    @contact: mathieu.doucet@nist.gov
 """
 
@@ -16,5 +28,5 @@
 class EllipticalCylinderModel(CEllipticalCylinderModel, BaseComponent):
     """ Class that evaluates a EllipticalCylinderModel model. 
-    	This file was auto-generated from elliptical_cylinder.h.
+    	This file was auto-generated from ../c_extensions/elliptical_cylinder.h.
     	Refer to that file and the structure it contains
     	for details of the model.
@@ -76,4 +88,13 @@
         
         return CEllipticalCylinderModel.runXY(self, x)
+        
+    def set_dispersion(self, parameter, dispersion):
+        """
+            Set the dispersion object for a model parameter
+            @param parameter: name of the parameter [string]
+            @dispersion: dispersion object of type DispersionModel
+        """
+        return CEllipticalCylinderModel.set_dispersion(self, parameter, dispersion.cdisp)
+        
    
 # End of file
Index: sansmodels/src/sans/models/SphereModel.py
===================================================================
--- sansmodels/src/sans/models/SphereModel.py	(revision ae3ce4eaa593202ebf49122e196897c0d45a25b6)
+++ sansmodels/src/sans/models/SphereModel.py	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
@@ -1,11 +1,23 @@
 #!/usr/bin/env python
+"""
+	This software was developed by the University of Tennessee as part of the
+	Distributed Data Analysis of Neutron Scattering Experiments (DANSE)
+	project funded by the US National Science Foundation.
+
+	If you use DANSE applications to do scientific research that leads to
+	publication, we ask that you acknowledge the use of the software with the
+	following sentence:
+
+	"This work benefited from DANSE software developed under NSF award DMR-0520547."
+
+	copyright 2008, University of Tennessee
+"""
+
 """ Provide functionality for a C extension model
 
 	WARNING: THIS FILE WAS GENERATED BY WRAPPERGENERATOR.PY
- 	         DO NOT MODIFY THIS FILE, MODIFY sphere.h
+ 	         DO NOT MODIFY THIS FILE, MODIFY ../c_extensions/sphere.h
  	         AND RE-RUN THE GENERATOR SCRIPT
 
-    @author: Mathieu Doucet / UTK
-    @contact: mathieu.doucet@nist.gov
 """
 
@@ -16,5 +28,5 @@
 class SphereModel(CSphereModel, BaseComponent):
     """ Class that evaluates a SphereModel model. 
-    	This file was auto-generated from sphere.h.
+    	This file was auto-generated from ../c_extensions/sphere.h.
     	Refer to that file and the structure it contains
     	for details of the model.
@@ -66,4 +78,13 @@
         
         return CSphereModel.runXY(self, x)
+        
+    def set_dispersion(self, parameter, dispersion):
+        """
+            Set the dispersion object for a model parameter
+            @param parameter: name of the parameter [string]
+            @dispersion: dispersion object of type DispersionModel
+        """
+        return CSphereModel.set_dispersion(self, parameter, dispersion.cdisp)
+        
    
 # End of file
Index: sansmodels/src/sans/models/c_extensions/c_models.cpp
===================================================================
--- sansmodels/src/sans/models/c_extensions/c_models.cpp	(revision af03ddd0932915d6659404ab7824e4223ea2944f)
+++ sansmodels/src/sans/models/c_extensions/c_models.cpp	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
@@ -9,10 +9,16 @@
 
 void addCCylinderModel(PyObject *module);
+void addCSphereModel(PyObject *module);
+void addCCoreShellModel(PyObject *module);
+void addCCoreShellCylinderModel(PyObject *module);
+void addCEllipsoidModel(PyObject *module);
+void addCEllipticalCylinderModel(PyObject *module);
+
 extern "C" {
-	void addCCoreShellCylinderModel(PyObject *module);
-	void addCCoreShellModel(PyObject *module);
-	void addCEllipsoidModel(PyObject *module);
-	void addCSphereModel(PyObject *module);
-	void addCEllipticalCylinderModel(PyObject *module);
+	//void addCCoreShellCylinderModel(PyObject *module);
+	//void addCCoreShellModel(PyObject *module);
+	//void addCEllipsoidModel(PyObject *module);
+	//void addCSphereModel(PyObject *module);
+	//void addCEllipticalCylinderModel(PyObject *module);
 	void addDisperser(PyObject *module);
 	void addCGaussian(PyObject *module);
Index: sansmodels/src/sans/models/c_extensions/core_shell.h
===================================================================
--- sansmodels/src/sans/models/c_extensions/core_shell.h	(revision ae3ce4eaa593202ebf49122e196897c0d45a25b6)
+++ sansmodels/src/sans/models/c_extensions/core_shell.h	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
@@ -2,10 +2,11 @@
 #define core_shell_h
 
-/** 
- * Structure definition for core-shell parameters 
+/**
+ * Structure definition for core-shell parameters
  */
  //[PYTHONCLASS] = CoreShellModel
+ //[DISP_PARAMS] = radius, thickness
 typedef struct {
-    /// Scale factor 
+    /// Scale factor
     //  [DEFAULT]=scale=1.0
     double scale;
@@ -27,5 +28,5 @@
 	/// Incoherent Background (cm-1) 0.000
 	//  [DEFAULT]=background=0 cm-1
-	double background;    
+	double background;
 } CoreShellParameters;
 
Index: sansmodels/src/sans/models/c_extensions/core_shell_cylinder.h
===================================================================
--- sansmodels/src/sans/models/c_extensions/core_shell_cylinder.h	(revision ae3ce4eaa593202ebf49122e196897c0d45a25b6)
+++ sansmodels/src/sans/models/c_extensions/core_shell_cylinder.h	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
@@ -2,49 +2,50 @@
 #define core_shell_cylinder_h
 
-/** 
- * Structure definition for core-shell cylinder parameters 
+/**
+ * Structure definition for core-shell cylinder parameters
  */
  //[PYTHONCLASS] = CoreShellCylinderModel
+ //[DISP_PARAMS] = radius, thickness, length, axis_theta, axis_phi
 typedef struct {
-    /// Scale factor 
+    /// Scale factor
     //  [DEFAULT]=scale=1.0
     double scale;
-    
+
     /// Core radius [A]
     //  [DEFAULT]=radius=20.0 A
     double radius;
-    
+
     /// Shell thickness [A]
     //  [DEFAULT]=thickness=10.0 A
     double thickness;
-    
+
     /// Core length [A]
     //  [DEFAULT]=length=400.0 A
     double length;
-    
+
     /// Core SLD [A-2]
     //  [DEFAULT]=core_sld=1.0e-6 A-2
     double core_sld;
-    
+
     /// Shell SLD [A-2]
     //  [DEFAULT]=shell_sld=4.0e-6 A-2
     double shell_sld;
-    
+
     /// Solvent SLD [A-2]
     //  [DEFAULT]=solvent_sld=1.0e-6 A-2
     double solvent_sld;
-    
+
 	/// Incoherent Background [cm-1]
 	//  [DEFAULT]=background=0 cm-1
-	double background;    
-	
+	double background;
+
     /// Orientation of the long axis of the core-shell cylinder w/respect incoming beam [rad]
     //  [DEFAULT]=axis_theta=1.57 rad
     double axis_theta;
-    
+
     /// Orientation of the long axis of the core-shell cylinder in the plane of the detector [rad]
     //  [DEFAULT]=axis_phi=0.0 rad
-    double axis_phi;    
-	
+    double axis_phi;
+
 } CoreShellCylinderParameters;
 
Index: sansmodels/src/sans/models/c_extensions/ellipsoid.h
===================================================================
--- sansmodels/src/sans/models/c_extensions/ellipsoid.h	(revision ae3ce4eaa593202ebf49122e196897c0d45a25b6)
+++ sansmodels/src/sans/models/c_extensions/ellipsoid.h	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
@@ -2,31 +2,32 @@
 #define ellipsoid_h
 
-/** 
- * Structure definition for ellipsoid parameters 
+/**
+ * Structure definition for ellipsoid parameters
  * The ellipsoid has axes radius_b, radius_b, radius_a.
  * Ref: Jan Skov Pedersen, Advances in Colloid and Interface Science, 70 (1997) 171-210
  */
  //[PYTHONCLASS] = EllipsoidModel
+ //[DISP_PARAMS] = radius_a, radius_b, axis_theta, axis_phi
 typedef struct {
-    /// Scale factor 
+    /// Scale factor
     //  [DEFAULT]=scale=1.0
     double scale;
-    
-    ///	Rotation axis radius_a [A] 
+
+    ///	Rotation axis radius_a [A]
     //  [DEFAULT]=radius_a=20.0 A
     double radius_a;
-    
-    /// Radius_b [A] 
+
+    /// Radius_b [A]
     //  [DEFAULT]=radius_b=400 A
     double radius_b;
-    
-    ///	Contrast [Å-2] 
+
+    ///	Contrast [Å-2]
     //  [DEFAULT]=contrast=3.0e-6 A-2
     double contrast;
-    
+
 	/// Incoherent Background [cm-1]
 	//  [DEFAULT]=background=0 cm-1
-	double background;    
-	
+	double background;
+
     /// Orientation of the long axis of the ellipsoid w/respect incoming beam [rad]
     //  [DEFAULT]=axis_theta=1.57 rad
@@ -34,5 +35,5 @@
     /// Orientation of the long axis of the ellipsoid in the plane of the detector [rad]
     //  [DEFAULT]=axis_phi=0.0 rad
-    double axis_phi;    
+    double axis_phi;
 } EllipsoidParameters;
 
Index: sansmodels/src/sans/models/c_extensions/elliptical_cylinder.h
===================================================================
--- sansmodels/src/sans/models/c_extensions/elliptical_cylinder.h	(revision ae3ce4eaa593202ebf49122e196897c0d45a25b6)
+++ sansmodels/src/sans/models/c_extensions/elliptical_cylinder.h	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
@@ -2,9 +2,11 @@
 #define ell_cylinder_h
 
-/** Structure definition for cylinder parameters 
+/** Structure definition for cylinder parameters
  * [PYTHONCLASS] = EllipticalCylinderModel
+ * [DISP_PARAMS] = r_minor, r_ratio, length, cyl_theta, cyl_phi, cyl_psi
+ *
  * */
 typedef struct {
-    /// Scale factor 
+    /// Scale factor
     //  [DEFAULT]=scale=1.0
     double scale;
@@ -23,5 +25,5 @@
 	/// Incoherent Background (cm-1) 0.000
 	//  [DEFAULT]=background=0 cm-1
-	double background;    
+	double background;
     /// Orientation of the cylinder axis w/respect incoming beam [rad]
     //  [DEFAULT]=cyl_theta=1.57 rad
@@ -29,8 +31,8 @@
     /// Orientation of the cylinder in the plane of the detector [rad]
     //  [DEFAULT]=cyl_phi=0.0 rad
-    double cyl_phi;    
+    double cyl_phi;
     /// Orientation of major radius of the cross-section w/respect vector q [rad]
     //  [DEFAULT]=cyl_psi=0.0 rad
-    double cyl_psi;    
+    double cyl_psi;
 } EllipticalCylinderParameters;
 
Index: sansmodels/src/sans/models/c_extensions/sphere.h
===================================================================
--- sansmodels/src/sans/models/c_extensions/sphere.h	(revision ae3ce4eaa593202ebf49122e196897c0d45a25b6)
+++ sansmodels/src/sans/models/c_extensions/sphere.h	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
@@ -2,24 +2,25 @@
 #define sphere_h
 
-/** 
- * Structure definition for sphere parameters 
+/**
+ * Structure definition for sphere parameters
  */
  //[PYTHONCLASS] = SphereModel
+ //[DISP_PARAMS] = radius
 typedef struct {
-    /// Scale factor 
+    /// Scale factor
     //  [DEFAULT]=scale=1.0e-6
     double scale;
-    
-    ///	Radius of sphere [A] 
+
+    ///	Radius of sphere [A]
     //  [DEFAULT]=radius=60.0 A
     double radius;
-    
-    ///	Contrast [Å-2] 
+
+    ///	Contrast [Å-2]
     //  [DEFAULT]=contrast=1.0 A-2
     double contrast;
-    
+
 	/// Incoherent Background [cm-1]
 	//  [DEFAULT]=background=0 cm-1
-	double background;    
+	double background;
 } SphereParameters;
 
Index: sansmodels/src/sans/models/c_models/CCoreShellCylinderModel.cpp
===================================================================
--- sansmodels/src/sans/models/c_models/CCoreShellCylinderModel.cpp	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
+++ sansmodels/src/sans/models/c_models/CCoreShellCylinderModel.cpp	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
@@ -0,0 +1,417 @@
+/**
+	This software was developed by the University of Tennessee as part of the
+	Distributed Data Analysis of Neutron Scattering Experiments (DANSE)
+	project funded by the US National Science Foundation.
+
+	If you use DANSE applications to do scientific research that leads to
+	publication, we ask that you acknowledge the use of the software with the
+	following sentence:
+
+	"This work benefited from DANSE software developed under NSF award DMR-0520547."
+
+	copyright 2008, University of Tennessee
+ */
+
+/** CCoreShellCylinderModel
+ *
+ * C extension 
+ *
+ * WARNING: THIS FILE WAS GENERATED BY WRAPPERGENERATOR.PY
+ *          DO NOT MODIFY THIS FILE, MODIFY core_shell_cylinder.h
+ *          AND RE-RUN THE GENERATOR SCRIPT
+ *
+ */
+ 
+extern "C" {
+#include <Python.h>
+#include "structmember.h"
+#include <stdio.h>
+#include <stdlib.h>
+#include <math.h>
+#include <time.h>
+#include "core_shell_cylinder.h"
+}
+
+#include "models.hh"
+#include "dispersion_visitor.hh"
+
+/// Error object for raised exceptions
+static PyObject * CCoreShellCylinderModelError = NULL;
+
+
+// Class definition
+typedef struct {
+    PyObject_HEAD
+    /// Parameters
+    PyObject * params;
+    /// Dispersion parameters
+    PyObject * dispersion;
+    /// Underlying model object
+    CoreShellCylinderModel * model;
+    /// Log for unit testing
+    PyObject * log;
+} CCoreShellCylinderModel;
+
+
+static void
+CCoreShellCylinderModel_dealloc(CCoreShellCylinderModel* self)
+{
+    self->ob_type->tp_free((PyObject*)self);
+    
+
+}
+
+static PyObject *
+CCoreShellCylinderModel_new(PyTypeObject *type, PyObject *args, PyObject *kwds)
+{
+    CCoreShellCylinderModel *self;
+    
+    self = (CCoreShellCylinderModel *)type->tp_alloc(type, 0);
+   
+    return (PyObject *)self;
+}
+
+static int
+CCoreShellCylinderModel_init(CCoreShellCylinderModel *self, PyObject *args, PyObject *kwds)
+{
+    if (self != NULL) {
+    	
+    	// Create parameters
+        self->params = PyDict_New();
+        self->dispersion = PyDict_New();
+        self->model = new CoreShellCylinderModel();
+        
+        // Initialize parameter dictionary
+        PyDict_SetItemString(self->params,"scale",Py_BuildValue("d",1.000000));
+        PyDict_SetItemString(self->params,"axis_theta",Py_BuildValue("d",1.570000));
+        PyDict_SetItemString(self->params,"solvent_sld",Py_BuildValue("d",0.000001));
+        PyDict_SetItemString(self->params,"thickness",Py_BuildValue("d",10.000000));
+        PyDict_SetItemString(self->params,"axis_phi",Py_BuildValue("d",0.000000));
+        PyDict_SetItemString(self->params,"length",Py_BuildValue("d",400.000000));
+        PyDict_SetItemString(self->params,"core_sld",Py_BuildValue("d",0.000001));
+        PyDict_SetItemString(self->params,"radius",Py_BuildValue("d",20.000000));
+        PyDict_SetItemString(self->params,"background",Py_BuildValue("d",0.000000));
+        PyDict_SetItemString(self->params,"shell_sld",Py_BuildValue("d",0.000004));
+        // Initialize dispersion / averaging parameter dict
+        DispersionVisitor* visitor = new DispersionVisitor();
+        PyObject * disp_dict;
+        disp_dict = PyDict_New();
+        self->model->radius.dispersion->accept_as_source(visitor, self->model->radius.dispersion, disp_dict);
+        PyDict_SetItemString(self->dispersion, "radius", disp_dict);
+        disp_dict = PyDict_New();
+        self->model->thickness.dispersion->accept_as_source(visitor, self->model->thickness.dispersion, disp_dict);
+        PyDict_SetItemString(self->dispersion, "thickness", disp_dict);
+        disp_dict = PyDict_New();
+        self->model->length.dispersion->accept_as_source(visitor, self->model->length.dispersion, disp_dict);
+        PyDict_SetItemString(self->dispersion, "length", disp_dict);
+        disp_dict = PyDict_New();
+        self->model->axis_theta.dispersion->accept_as_source(visitor, self->model->axis_theta.dispersion, disp_dict);
+        PyDict_SetItemString(self->dispersion, "axis_theta", disp_dict);
+        disp_dict = PyDict_New();
+        self->model->axis_phi.dispersion->accept_as_source(visitor, self->model->axis_phi.dispersion, disp_dict);
+        PyDict_SetItemString(self->dispersion, "axis_phi", disp_dict);
+
+
+         
+        // Create empty log
+        self->log = PyDict_New();
+        
+        
+
+    }
+    return 0;
+}
+
+static PyMemberDef CCoreShellCylinderModel_members[] = {
+    {"params", T_OBJECT, offsetof(CCoreShellCylinderModel, params), 0,
+     "Parameters"},
+	{"dispersion", T_OBJECT, offsetof(CCoreShellCylinderModel, dispersion), 0,
+	  "Dispersion parameters"},     
+    {"log", T_OBJECT, offsetof(CCoreShellCylinderModel, log), 0,
+     "Log"},
+    {NULL}  /* Sentinel */
+};
+
+/** Read double from PyObject
+    @param p PyObject
+    @return double
+*/
+double CCoreShellCylinderModel_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(CCoreShellCylinderModel *self, PyObject *args) {
+	double q_value, phi_value;
+	PyObject* pars;
+	int npars;
+	
+	// Get parameters
+	
+	    // Reader parameter dictionary
+    self->model->scale = PyFloat_AsDouble( PyDict_GetItemString(self->params, "scale") );
+    self->model->axis_theta = PyFloat_AsDouble( PyDict_GetItemString(self->params, "axis_theta") );
+    self->model->solvent_sld = PyFloat_AsDouble( PyDict_GetItemString(self->params, "solvent_sld") );
+    self->model->thickness = PyFloat_AsDouble( PyDict_GetItemString(self->params, "thickness") );
+    self->model->axis_phi = PyFloat_AsDouble( PyDict_GetItemString(self->params, "axis_phi") );
+    self->model->length = PyFloat_AsDouble( PyDict_GetItemString(self->params, "length") );
+    self->model->core_sld = PyFloat_AsDouble( PyDict_GetItemString(self->params, "core_sld") );
+    self->model->radius = PyFloat_AsDouble( PyDict_GetItemString(self->params, "radius") );
+    self->model->background = PyFloat_AsDouble( PyDict_GetItemString(self->params, "background") );
+    self->model->shell_sld = PyFloat_AsDouble( PyDict_GetItemString(self->params, "shell_sld") );
+    // Read in dispersion parameters
+    PyObject* disp_dict;
+    DispersionVisitor* visitor = new DispersionVisitor();
+    disp_dict = PyDict_GetItemString(self->dispersion, "radius");
+    self->model->radius.dispersion->accept_as_destination(visitor, self->model->radius.dispersion, disp_dict);
+    disp_dict = PyDict_GetItemString(self->dispersion, "thickness");
+    self->model->thickness.dispersion->accept_as_destination(visitor, self->model->thickness.dispersion, disp_dict);
+    disp_dict = PyDict_GetItemString(self->dispersion, "length");
+    self->model->length.dispersion->accept_as_destination(visitor, self->model->length.dispersion, disp_dict);
+    disp_dict = PyDict_GetItemString(self->dispersion, "axis_theta");
+    self->model->axis_theta.dispersion->accept_as_destination(visitor, self->model->axis_theta.dispersion, disp_dict);
+    disp_dict = PyDict_GetItemString(self->dispersion, "axis_phi");
+    self->model->axis_phi.dispersion->accept_as_destination(visitor, self->model->axis_phi.dispersion, disp_dict);
+
+	
+	// Get input and determine whether we have to supply a 1D or 2D return value.
+	if ( !PyArg_ParseTuple(args,"O",&pars) ) {
+	    PyErr_SetString(CCoreShellCylinderModelError, 
+	    	"CCoreShellCylinderModel.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(CCoreShellCylinderModelError, 
+	    		"CCoreShellCylinderModel.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 = CCoreShellCylinderModel_readDouble(PyList_GET_ITEM(pars,0));
+	    phi_value = CCoreShellCylinderModel_readDouble(PyList_GET_ITEM(pars,1));
+	    // Skip zero
+	    if (q_value==0) {
+	    	return Py_BuildValue("d",0.0);
+	    }
+		return Py_BuildValue("d",(*(self->model)).evaluate_rphi(q_value,phi_value));
+
+	} else {
+
+		// We have a scalar q, we will evaluate I(q)
+		q_value = CCoreShellCylinderModel_readDouble(pars);		
+		
+		return Py_BuildValue("d",(*(self->model))(q_value));
+	}	
+}
+
+/**
+ * Function to call to evaluate model in cartesian coordinates
+ * @param args: input q or [qx, qy]]
+ * @return: function value
+ */
+static PyObject * runXY(CCoreShellCylinderModel *self, PyObject *args) {
+	double qx_value, qy_value;
+	PyObject* pars;
+	int npars;
+	
+	// Get parameters
+	
+	    // Reader parameter dictionary
+    self->model->scale = PyFloat_AsDouble( PyDict_GetItemString(self->params, "scale") );
+    self->model->axis_theta = PyFloat_AsDouble( PyDict_GetItemString(self->params, "axis_theta") );
+    self->model->solvent_sld = PyFloat_AsDouble( PyDict_GetItemString(self->params, "solvent_sld") );
+    self->model->thickness = PyFloat_AsDouble( PyDict_GetItemString(self->params, "thickness") );
+    self->model->axis_phi = PyFloat_AsDouble( PyDict_GetItemString(self->params, "axis_phi") );
+    self->model->length = PyFloat_AsDouble( PyDict_GetItemString(self->params, "length") );
+    self->model->core_sld = PyFloat_AsDouble( PyDict_GetItemString(self->params, "core_sld") );
+    self->model->radius = PyFloat_AsDouble( PyDict_GetItemString(self->params, "radius") );
+    self->model->background = PyFloat_AsDouble( PyDict_GetItemString(self->params, "background") );
+    self->model->shell_sld = PyFloat_AsDouble( PyDict_GetItemString(self->params, "shell_sld") );
+    // Read in dispersion parameters
+    PyObject* disp_dict;
+    DispersionVisitor* visitor = new DispersionVisitor();
+    disp_dict = PyDict_GetItemString(self->dispersion, "radius");
+    self->model->radius.dispersion->accept_as_destination(visitor, self->model->radius.dispersion, disp_dict);
+    disp_dict = PyDict_GetItemString(self->dispersion, "thickness");
+    self->model->thickness.dispersion->accept_as_destination(visitor, self->model->thickness.dispersion, disp_dict);
+    disp_dict = PyDict_GetItemString(self->dispersion, "length");
+    self->model->length.dispersion->accept_as_destination(visitor, self->model->length.dispersion, disp_dict);
+    disp_dict = PyDict_GetItemString(self->dispersion, "axis_theta");
+    self->model->axis_theta.dispersion->accept_as_destination(visitor, self->model->axis_theta.dispersion, disp_dict);
+    disp_dict = PyDict_GetItemString(self->dispersion, "axis_phi");
+    self->model->axis_phi.dispersion->accept_as_destination(visitor, self->model->axis_phi.dispersion, disp_dict);
+
+	
+	// Get input and determine whether we have to supply a 1D or 2D return value.
+	if ( !PyArg_ParseTuple(args,"O",&pars) ) {
+	    PyErr_SetString(CCoreShellCylinderModelError, 
+	    	"CCoreShellCylinderModel.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(CCoreShellCylinderModelError, 
+	    		"CCoreShellCylinderModel.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 = CCoreShellCylinderModel_readDouble(PyList_GET_ITEM(pars,0));
+	    qy_value = CCoreShellCylinderModel_readDouble(PyList_GET_ITEM(pars,1));
+	    return Py_BuildValue("d",(*(self->model))(qx_value,qy_value));
+
+	} else {
+
+		// We have a scalar q, we will evaluate I(q)
+		qx_value = CCoreShellCylinderModel_readDouble(pars);		
+		
+		return Py_BuildValue("d",(*(self->model))(qx_value));
+	}	
+}
+
+static PyObject * reset(CCoreShellCylinderModel *self, PyObject *args) {
+    
+
+    return Py_BuildValue("d",0.0);
+}
+
+static PyObject * set_dispersion(CCoreShellCylinderModel *self, PyObject *args) {
+	PyObject * disp;
+	const char * par_name;
+
+	if ( !PyArg_ParseTuple(args,"sO", &par_name, &disp) ) {
+	    PyErr_SetString(CCoreShellCylinderModelError,
+	    	"CCoreShellCylinderModel.set_dispersion expects a DispersionModel object.");
+		return NULL;
+	}
+	void *temp = PyCObject_AsVoidPtr(disp);
+	DispersionModel * dispersion = static_cast<DispersionModel *>(temp);
+
+
+	// Ugliness necessary to go from python to C
+	    // TODO: refactor this
+    if (!strcmp(par_name, "radius")) {
+        self->model->radius.dispersion = dispersion;
+    } else    if (!strcmp(par_name, "thickness")) {
+        self->model->thickness.dispersion = dispersion;
+    } else    if (!strcmp(par_name, "length")) {
+        self->model->length.dispersion = dispersion;
+    } else    if (!strcmp(par_name, "axis_theta")) {
+        self->model->axis_theta.dispersion = dispersion;
+    } else    if (!strcmp(par_name, "axis_phi")) {
+        self->model->axis_phi.dispersion = dispersion;
+    } else {
+	    PyErr_SetString(CCoreShellCylinderModelError,
+	    	"CCoreShellCylinderModel.set_dispersion expects a valid parameter name.");
+		return NULL;
+	}
+
+	DispersionVisitor* visitor = new DispersionVisitor();
+	PyObject * disp_dict = PyDict_New();
+	dispersion->accept_as_source(visitor, dispersion, disp_dict);
+	PyDict_SetItemString(self->dispersion, par_name, disp_dict);
+    return Py_BuildValue("i",1);
+}
+
+
+static PyMethodDef CCoreShellCylinderModel_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"},
+    {"set_dispersion",      (PyCFunction)set_dispersion     , METH_VARARGS,
+      "Set the dispersion model for a given parameter"},
+   {NULL}
+};
+
+static PyTypeObject CCoreShellCylinderModelType = {
+    PyObject_HEAD_INIT(NULL)
+    0,                         /*ob_size*/
+    "CCoreShellCylinderModel",             /*tp_name*/
+    sizeof(CCoreShellCylinderModel),             /*tp_basicsize*/
+    0,                         /*tp_itemsize*/
+    (destructor)CCoreShellCylinderModel_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*/
+    "CCoreShellCylinderModel objects",           /* tp_doc */
+    0,		               /* tp_traverse */
+    0,		               /* tp_clear */
+    0,		               /* tp_richcompare */
+    0,		               /* tp_weaklistoffset */
+    0,		               /* tp_iter */
+    0,		               /* tp_iternext */
+    CCoreShellCylinderModel_methods,             /* tp_methods */
+    CCoreShellCylinderModel_members,             /* tp_members */
+    0,                         /* tp_getset */
+    0,                         /* tp_base */
+    0,                         /* tp_dict */
+    0,                         /* tp_descr_get */
+    0,                         /* tp_descr_set */
+    0,                         /* tp_dictoffset */
+    (initproc)CCoreShellCylinderModel_init,      /* tp_init */
+    0,                         /* tp_alloc */
+    CCoreShellCylinderModel_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 addCCoreShellCylinderModel(PyObject *module) {
+	PyObject *d;
+	
+    if (PyType_Ready(&CCoreShellCylinderModelType) < 0)
+        return;
+
+    Py_INCREF(&CCoreShellCylinderModelType);
+    PyModule_AddObject(module, "CCoreShellCylinderModel", (PyObject *)&CCoreShellCylinderModelType);
+    
+    d = PyModule_GetDict(module);
+    CCoreShellCylinderModelError = PyErr_NewException("CCoreShellCylinderModel.error", NULL, NULL);
+    PyDict_SetItemString(d, "CCoreShellCylinderModelError", CCoreShellCylinderModelError);
+}
+
Index: sansmodels/src/sans/models/c_models/CCoreShellModel.cpp
===================================================================
--- sansmodels/src/sans/models/c_models/CCoreShellModel.cpp	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
+++ sansmodels/src/sans/models/c_models/CCoreShellModel.cpp	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
@@ -0,0 +1,381 @@
+/**
+	This software was developed by the University of Tennessee as part of the
+	Distributed Data Analysis of Neutron Scattering Experiments (DANSE)
+	project funded by the US National Science Foundation.
+
+	If you use DANSE applications to do scientific research that leads to
+	publication, we ask that you acknowledge the use of the software with the
+	following sentence:
+
+	"This work benefited from DANSE software developed under NSF award DMR-0520547."
+
+	copyright 2008, University of Tennessee
+ */
+
+/** CCoreShellModel
+ *
+ * C extension 
+ *
+ * WARNING: THIS FILE WAS GENERATED BY WRAPPERGENERATOR.PY
+ *          DO NOT MODIFY THIS FILE, MODIFY core_shell.h
+ *          AND RE-RUN THE GENERATOR SCRIPT
+ *
+ */
+ 
+extern "C" {
+#include <Python.h>
+#include "structmember.h"
+#include <stdio.h>
+#include <stdlib.h>
+#include <math.h>
+#include <time.h>
+#include "core_shell.h"
+}
+
+#include "models.hh"
+#include "dispersion_visitor.hh"
+
+/// Error object for raised exceptions
+static PyObject * CCoreShellModelError = NULL;
+
+
+// Class definition
+typedef struct {
+    PyObject_HEAD
+    /// Parameters
+    PyObject * params;
+    /// Dispersion parameters
+    PyObject * dispersion;
+    /// Underlying model object
+    CoreShellModel * model;
+    /// Log for unit testing
+    PyObject * log;
+} CCoreShellModel;
+
+
+static void
+CCoreShellModel_dealloc(CCoreShellModel* self)
+{
+    self->ob_type->tp_free((PyObject*)self);
+    
+
+}
+
+static PyObject *
+CCoreShellModel_new(PyTypeObject *type, PyObject *args, PyObject *kwds)
+{
+    CCoreShellModel *self;
+    
+    self = (CCoreShellModel *)type->tp_alloc(type, 0);
+   
+    return (PyObject *)self;
+}
+
+static int
+CCoreShellModel_init(CCoreShellModel *self, PyObject *args, PyObject *kwds)
+{
+    if (self != NULL) {
+    	
+    	// Create parameters
+        self->params = PyDict_New();
+        self->dispersion = PyDict_New();
+        self->model = new CoreShellModel();
+        
+        // Initialize parameter dictionary
+        PyDict_SetItemString(self->params,"core_sld",Py_BuildValue("d",0.000001));
+        PyDict_SetItemString(self->params,"thickness",Py_BuildValue("d",10.000000));
+        PyDict_SetItemString(self->params,"solvent_sld",Py_BuildValue("d",0.000003));
+        PyDict_SetItemString(self->params,"scale",Py_BuildValue("d",1.000000));
+        PyDict_SetItemString(self->params,"radius",Py_BuildValue("d",60.000000));
+        PyDict_SetItemString(self->params,"background",Py_BuildValue("d",0.000000));
+        PyDict_SetItemString(self->params,"shell_sld",Py_BuildValue("d",0.000002));
+        // Initialize dispersion / averaging parameter dict
+        DispersionVisitor* visitor = new DispersionVisitor();
+        PyObject * disp_dict;
+        disp_dict = PyDict_New();
+        self->model->radius.dispersion->accept_as_source(visitor, self->model->radius.dispersion, disp_dict);
+        PyDict_SetItemString(self->dispersion, "radius", disp_dict);
+        disp_dict = PyDict_New();
+        self->model->thickness.dispersion->accept_as_source(visitor, self->model->thickness.dispersion, disp_dict);
+        PyDict_SetItemString(self->dispersion, "thickness", disp_dict);
+
+
+         
+        // Create empty log
+        self->log = PyDict_New();
+        
+        
+
+    }
+    return 0;
+}
+
+static PyMemberDef CCoreShellModel_members[] = {
+    {"params", T_OBJECT, offsetof(CCoreShellModel, params), 0,
+     "Parameters"},
+	{"dispersion", T_OBJECT, offsetof(CCoreShellModel, dispersion), 0,
+	  "Dispersion parameters"},     
+    {"log", T_OBJECT, offsetof(CCoreShellModel, log), 0,
+     "Log"},
+    {NULL}  /* Sentinel */
+};
+
+/** Read double from PyObject
+    @param p PyObject
+    @return double
+*/
+double CCoreShellModel_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(CCoreShellModel *self, PyObject *args) {
+	double q_value, phi_value;
+	PyObject* pars;
+	int npars;
+	
+	// Get parameters
+	
+	    // Reader parameter dictionary
+    self->model->core_sld = PyFloat_AsDouble( PyDict_GetItemString(self->params, "core_sld") );
+    self->model->thickness = PyFloat_AsDouble( PyDict_GetItemString(self->params, "thickness") );
+    self->model->solvent_sld = PyFloat_AsDouble( PyDict_GetItemString(self->params, "solvent_sld") );
+    self->model->scale = PyFloat_AsDouble( PyDict_GetItemString(self->params, "scale") );
+    self->model->radius = PyFloat_AsDouble( PyDict_GetItemString(self->params, "radius") );
+    self->model->background = PyFloat_AsDouble( PyDict_GetItemString(self->params, "background") );
+    self->model->shell_sld = PyFloat_AsDouble( PyDict_GetItemString(self->params, "shell_sld") );
+    // Read in dispersion parameters
+    PyObject* disp_dict;
+    DispersionVisitor* visitor = new DispersionVisitor();
+    disp_dict = PyDict_GetItemString(self->dispersion, "radius");
+    self->model->radius.dispersion->accept_as_destination(visitor, self->model->radius.dispersion, disp_dict);
+    disp_dict = PyDict_GetItemString(self->dispersion, "thickness");
+    self->model->thickness.dispersion->accept_as_destination(visitor, self->model->thickness.dispersion, disp_dict);
+
+	
+	// Get input and determine whether we have to supply a 1D or 2D return value.
+	if ( !PyArg_ParseTuple(args,"O",&pars) ) {
+	    PyErr_SetString(CCoreShellModelError, 
+	    	"CCoreShellModel.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(CCoreShellModelError, 
+	    		"CCoreShellModel.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 = CCoreShellModel_readDouble(PyList_GET_ITEM(pars,0));
+	    phi_value = CCoreShellModel_readDouble(PyList_GET_ITEM(pars,1));
+	    // Skip zero
+	    if (q_value==0) {
+	    	return Py_BuildValue("d",0.0);
+	    }
+		return Py_BuildValue("d",(*(self->model)).evaluate_rphi(q_value,phi_value));
+
+	} else {
+
+		// We have a scalar q, we will evaluate I(q)
+		q_value = CCoreShellModel_readDouble(pars);		
+		
+		return Py_BuildValue("d",(*(self->model))(q_value));
+	}	
+}
+
+/**
+ * Function to call to evaluate model in cartesian coordinates
+ * @param args: input q or [qx, qy]]
+ * @return: function value
+ */
+static PyObject * runXY(CCoreShellModel *self, PyObject *args) {
+	double qx_value, qy_value;
+	PyObject* pars;
+	int npars;
+	
+	// Get parameters
+	
+	    // Reader parameter dictionary
+    self->model->core_sld = PyFloat_AsDouble( PyDict_GetItemString(self->params, "core_sld") );
+    self->model->thickness = PyFloat_AsDouble( PyDict_GetItemString(self->params, "thickness") );
+    self->model->solvent_sld = PyFloat_AsDouble( PyDict_GetItemString(self->params, "solvent_sld") );
+    self->model->scale = PyFloat_AsDouble( PyDict_GetItemString(self->params, "scale") );
+    self->model->radius = PyFloat_AsDouble( PyDict_GetItemString(self->params, "radius") );
+    self->model->background = PyFloat_AsDouble( PyDict_GetItemString(self->params, "background") );
+    self->model->shell_sld = PyFloat_AsDouble( PyDict_GetItemString(self->params, "shell_sld") );
+    // Read in dispersion parameters
+    PyObject* disp_dict;
+    DispersionVisitor* visitor = new DispersionVisitor();
+    disp_dict = PyDict_GetItemString(self->dispersion, "radius");
+    self->model->radius.dispersion->accept_as_destination(visitor, self->model->radius.dispersion, disp_dict);
+    disp_dict = PyDict_GetItemString(self->dispersion, "thickness");
+    self->model->thickness.dispersion->accept_as_destination(visitor, self->model->thickness.dispersion, disp_dict);
+
+	
+	// Get input and determine whether we have to supply a 1D or 2D return value.
+	if ( !PyArg_ParseTuple(args,"O",&pars) ) {
+	    PyErr_SetString(CCoreShellModelError, 
+	    	"CCoreShellModel.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(CCoreShellModelError, 
+	    		"CCoreShellModel.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 = CCoreShellModel_readDouble(PyList_GET_ITEM(pars,0));
+	    qy_value = CCoreShellModel_readDouble(PyList_GET_ITEM(pars,1));
+	    return Py_BuildValue("d",(*(self->model))(qx_value,qy_value));
+
+	} else {
+
+		// We have a scalar q, we will evaluate I(q)
+		qx_value = CCoreShellModel_readDouble(pars);		
+		
+		return Py_BuildValue("d",(*(self->model))(qx_value));
+	}	
+}
+
+static PyObject * reset(CCoreShellModel *self, PyObject *args) {
+    
+
+    return Py_BuildValue("d",0.0);
+}
+
+static PyObject * set_dispersion(CCoreShellModel *self, PyObject *args) {
+	PyObject * disp;
+	const char * par_name;
+
+	if ( !PyArg_ParseTuple(args,"sO", &par_name, &disp) ) {
+	    PyErr_SetString(CCoreShellModelError,
+	    	"CCoreShellModel.set_dispersion expects a DispersionModel object.");
+		return NULL;
+	}
+	void *temp = PyCObject_AsVoidPtr(disp);
+	DispersionModel * dispersion = static_cast<DispersionModel *>(temp);
+
+
+	// Ugliness necessary to go from python to C
+	    // TODO: refactor this
+    if (!strcmp(par_name, "radius")) {
+        self->model->radius.dispersion = dispersion;
+    } else    if (!strcmp(par_name, "thickness")) {
+        self->model->thickness.dispersion = dispersion;
+    } else {
+	    PyErr_SetString(CCoreShellModelError,
+	    	"CCoreShellModel.set_dispersion expects a valid parameter name.");
+		return NULL;
+	}
+
+	DispersionVisitor* visitor = new DispersionVisitor();
+	PyObject * disp_dict = PyDict_New();
+	dispersion->accept_as_source(visitor, dispersion, disp_dict);
+	PyDict_SetItemString(self->dispersion, par_name, disp_dict);
+    return Py_BuildValue("i",1);
+}
+
+
+static PyMethodDef CCoreShellModel_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"},
+    {"set_dispersion",      (PyCFunction)set_dispersion     , METH_VARARGS,
+      "Set the dispersion model for a given parameter"},
+   {NULL}
+};
+
+static PyTypeObject CCoreShellModelType = {
+    PyObject_HEAD_INIT(NULL)
+    0,                         /*ob_size*/
+    "CCoreShellModel",             /*tp_name*/
+    sizeof(CCoreShellModel),             /*tp_basicsize*/
+    0,                         /*tp_itemsize*/
+    (destructor)CCoreShellModel_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*/
+    "CCoreShellModel objects",           /* tp_doc */
+    0,		               /* tp_traverse */
+    0,		               /* tp_clear */
+    0,		               /* tp_richcompare */
+    0,		               /* tp_weaklistoffset */
+    0,		               /* tp_iter */
+    0,		               /* tp_iternext */
+    CCoreShellModel_methods,             /* tp_methods */
+    CCoreShellModel_members,             /* tp_members */
+    0,                         /* tp_getset */
+    0,                         /* tp_base */
+    0,                         /* tp_dict */
+    0,                         /* tp_descr_get */
+    0,                         /* tp_descr_set */
+    0,                         /* tp_dictoffset */
+    (initproc)CCoreShellModel_init,      /* tp_init */
+    0,                         /* tp_alloc */
+    CCoreShellModel_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 addCCoreShellModel(PyObject *module) {
+	PyObject *d;
+	
+    if (PyType_Ready(&CCoreShellModelType) < 0)
+        return;
+
+    Py_INCREF(&CCoreShellModelType);
+    PyModule_AddObject(module, "CCoreShellModel", (PyObject *)&CCoreShellModelType);
+    
+    d = PyModule_GetDict(module);
+    CCoreShellModelError = PyErr_NewException("CCoreShellModel.error", NULL, NULL);
+    PyDict_SetItemString(d, "CCoreShellModelError", CCoreShellModelError);
+}
+
Index: sansmodels/src/sans/models/c_models/CCylinderModel.cpp
===================================================================
--- sansmodels/src/sans/models/c_models/CCylinderModel.cpp	(revision af03ddd0932915d6659404ab7824e4223ea2944f)
+++ sansmodels/src/sans/models/c_models/CCylinderModel.cpp	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
@@ -83,11 +83,11 @@
         
         // Initialize parameter dictionary
-        PyDict_SetItemString(self->params,"scale",Py_BuildValue("d", self->model->scale()));
-        PyDict_SetItemString(self->params,"length",Py_BuildValue("d", self->model->length()));
-        PyDict_SetItemString(self->params,"cyl_theta",Py_BuildValue("d", self->model->cyl_theta()));
-        PyDict_SetItemString(self->params,"background",Py_BuildValue("d", self->model->background()));
-        PyDict_SetItemString(self->params,"radius",Py_BuildValue("d", self->model->radius()));
-        PyDict_SetItemString(self->params,"contrast",Py_BuildValue("d", self->model->contrast()));
-        PyDict_SetItemString(self->params,"cyl_phi",Py_BuildValue("d", self->model->cyl_phi()));
+        PyDict_SetItemString(self->params,"scale",Py_BuildValue("d",1.000000));
+        PyDict_SetItemString(self->params,"length",Py_BuildValue("d",400.000000));
+        PyDict_SetItemString(self->params,"cyl_theta",Py_BuildValue("d",1.000000));
+        PyDict_SetItemString(self->params,"background",Py_BuildValue("d",0.000000));
+        PyDict_SetItemString(self->params,"radius",Py_BuildValue("d",20.000000));
+        PyDict_SetItemString(self->params,"contrast",Py_BuildValue("d",0.000003));
+        PyDict_SetItemString(self->params,"cyl_phi",Py_BuildValue("d",1.000000));
         // Initialize dispersion / averaging parameter dict
         DispersionVisitor* visitor = new DispersionVisitor();
Index: sansmodels/src/sans/models/c_models/CCylinder_bck.cpp
===================================================================
--- sansmodels/src/sans/models/c_models/CCylinder_bck.cpp	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
+++ sansmodels/src/sans/models/c_models/CCylinder_bck.cpp	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
@@ -0,0 +1,413 @@
+/**
+	This software was developed by the University of Tennessee as part of the
+	Distributed Data Analysis of Neutron Scattering Experiments (DANSE)
+	project funded by the US National Science Foundation.
+
+	If you use DANSE applications to do scientific research that leads to
+	publication, we ask that you acknowledge the use of the software with the
+	following sentence:
+
+	"This work benefited from DANSE software developed under NSF award DMR-0520547."
+
+	copyright 2008, University of Tennessee
+ */
+
+/** CCylinderModel
+ *
+ * Base python class for cylinder model
+ *
+ */
+extern "C" {
+#include <Python.h>
+#include "structmember.h"
+#include <stdio.h>
+#include <stdlib.h>
+#include <math.h>
+#include <time.h>
+#include "cylinder.h"
+}
+
+
+#include "models.hh"
+#include "dispersion_visitor.hh"
+
+/// Error object for raised exceptions
+static PyObject * CCylinderModelError = NULL;
+
+
+// Class definition
+typedef struct {
+    PyObject_HEAD
+    /// Parameters
+    PyObject * params;
+    PyObject * dispersion;
+    CylinderModel * model;
+    /// Log for unit testing
+    PyObject * log;
+} CCylinderModel;
+
+
+static void
+CCylinderModel_dealloc(CCylinderModel* self)
+{
+    self->ob_type->tp_free((PyObject*)self);
+
+
+}
+
+static PyObject *
+CCylinderModel_new(PyTypeObject *type, PyObject *args, PyObject *kwds)
+{
+    CCylinderModel *self;
+
+    self = (CCylinderModel *)type->tp_alloc(type, 0);
+
+    return (PyObject *)self;
+}
+
+static int
+CCylinderModel_init(CCylinderModel *self, PyObject *args, PyObject *kwds)
+{
+    if (self != NULL) {
+    	printf("Hello cylinder\n");
+    	// Create parameters
+        self->params = PyDict_New();
+        self->dispersion = PyDict_New();
+        self->model = new CylinderModel();
+
+        // Initialize parameter dictionary
+        PyDict_SetItemString(self->params,"scale",     Py_BuildValue("d", self->model->scale()));
+        PyDict_SetItemString(self->params,"length",    Py_BuildValue("d", self->model->length()));
+        PyDict_SetItemString(self->params,"cyl_theta", Py_BuildValue("d", self->model->cyl_theta()));
+        PyDict_SetItemString(self->params,"background",Py_BuildValue("d", self->model->background()));
+        PyDict_SetItemString(self->params,"radius",    Py_BuildValue("d", self->model->radius()));
+        PyDict_SetItemString(self->params,"contrast",  Py_BuildValue("d", self->model->contrast()));
+        PyDict_SetItemString(self->params,"cyl_phi",   Py_BuildValue("d", self->model->cyl_phi()));
+
+        // Initialize dispersion / averaging parameter dict
+        DispersionVisitor* visitor = new DispersionVisitor();
+        PyObject * disp_dict = PyDict_New();
+        self->model->radius.dispersion->accept_as_source(visitor, self->model->radius.dispersion, disp_dict);
+        PyDict_SetItemString(self->dispersion, "radius", disp_dict);
+
+        disp_dict = PyDict_New();
+        self->model->length.dispersion->accept_as_source(visitor, self->model->length.dispersion, disp_dict);
+        PyDict_SetItemString(self->dispersion, "length", disp_dict);
+
+        disp_dict = PyDict_New();
+        self->model->cyl_phi.dispersion->accept_as_source(visitor, self->model->cyl_phi.dispersion, disp_dict);
+        PyDict_SetItemString(self->dispersion, "cyl_phi", disp_dict);
+
+        disp_dict = PyDict_New();
+        self->model->cyl_theta.dispersion->accept_as_source(visitor, self->model->cyl_theta.dispersion, disp_dict);
+        PyDict_SetItemString(self->dispersion, "cyl_theta", disp_dict);
+
+        // Create empty log
+        self->log = PyDict_New();
+
+    }
+    return 0;
+}
+
+static PyMemberDef CCylinderModel_members[] = {
+	{"params", T_OBJECT, offsetof(CCylinderModel, params), 0,
+	 "Parameters"},
+	{"dispersion", T_OBJECT, offsetof(CCylinderModel, dispersion), 0,
+	  "Dispersion parameters"},
+    {"log", T_OBJECT, offsetof(CCylinderModel, log), 0,
+     "Log"},
+    {NULL}  /* Sentinel */
+};
+
+/** Read double from PyObject
+    @param p PyObject
+    @return double
+*/
+double CCylinderModel_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;
+    }
+}
+
+/**
+ * Method to set the dispersion model for a parameter.
+ * We need to update the dispersion object of the parameter
+ * on the C++ side and update the dispersion dictionary on the python side.
+ *
+ * TODO: define read-only integers that will be data members of this
+ * class and allow the user to specify which model is needed with a code name.
+ *
+ * This method should take in a code name (FLAT, GAUSSIAN) and a string
+ * that represent the parameter to attach the dispersion model to.
+ *
+ */
+static PyObject * set_dispersion_model(CCylinderModel *self, PyObject *args) {
+
+}
+
+/**
+ * Function to call to evaluate model
+ * @param args: input q or [q,phi]
+ * @return: function value
+ */
+static PyObject * run(CCylinderModel *self, PyObject *args) {
+	double q_value, phi_value;
+	PyObject* pars;
+	int npars;
+
+	// Get parameters
+
+	// Reader parameter dictionary
+    self->model->scale = PyFloat_AsDouble( PyDict_GetItemString(self->params, "scale") );
+    self->model->length = PyFloat_AsDouble( PyDict_GetItemString(self->params, "length") );
+    self->model->cyl_theta = PyFloat_AsDouble( PyDict_GetItemString(self->params, "cyl_theta") );
+    self->model->background = PyFloat_AsDouble( PyDict_GetItemString(self->params, "background") );
+    self->model->radius = PyFloat_AsDouble( PyDict_GetItemString(self->params, "radius") );
+    self->model->contrast = PyFloat_AsDouble( PyDict_GetItemString(self->params, "contrast") );
+    self->model->cyl_phi = PyFloat_AsDouble( PyDict_GetItemString(self->params, "cyl_phi") );
+
+    // Read in dispersion parameters
+    PyObject* disp_dict;
+
+    DispersionVisitor* visitor = new DispersionVisitor();
+    disp_dict = PyDict_GetItemString(self->dispersion, "radius");
+    self->model->radius.dispersion->accept_as_destination(visitor, self->model->radius.dispersion, disp_dict);
+
+    disp_dict = PyDict_GetItemString(self->dispersion, "length");
+    self->model->length.dispersion->accept_as_destination(visitor, self->model->length.dispersion, disp_dict);
+
+    disp_dict = PyDict_GetItemString(self->dispersion, "cyl_theta");
+    self->model->cyl_theta.dispersion->accept_as_destination(visitor, self->model->cyl_theta.dispersion, disp_dict);
+
+    disp_dict = PyDict_GetItemString(self->dispersion, "cyl_phi");
+    self->model->cyl_phi.dispersion->accept_as_destination(visitor, self->model->cyl_phi.dispersion, disp_dict);
+
+
+	// Get input and determine whether we have to supply a 1D or 2D return value.
+	if ( !PyArg_ParseTuple(args,"O",&pars) ) {
+	    PyErr_SetString(CCylinderModelError,
+	    	"CCylinderModel.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(CCylinderModelError,
+	    		"CCylinderModel.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 = CCylinderModel_readDouble(PyList_GET_ITEM(pars,0));
+	    phi_value = CCylinderModel_readDouble(PyList_GET_ITEM(pars,1));
+	    // Skip zero
+	    if (q_value==0) {
+	    	return Py_BuildValue("d",0.0);
+	    }
+		return Py_BuildValue("d",(*(self->model)).evaluate_rphi(q_value,phi_value));
+
+	} else {
+
+		// We have a scalar q, we will evaluate I(q)
+		q_value = CCylinderModel_readDouble(pars);
+
+		return Py_BuildValue("d",(*(self->model))(q_value));
+	}
+}
+
+/**
+ * Function to call to evaluate model in cartesian coordinates
+ * @param args: input q or [qx, qy]]
+ * @return: function value
+ */
+static PyObject * runXY(CCylinderModel *self, PyObject *args) {
+	double qx_value, qy_value;
+	PyObject* pars;
+	int npars;
+
+	// Get parameters
+
+	// Reader parameter dictionary
+    self->model->scale = PyFloat_AsDouble( PyDict_GetItemString(self->params, "scale") );
+    self->model->length = PyFloat_AsDouble( PyDict_GetItemString(self->params, "length") );
+    self->model->cyl_theta = PyFloat_AsDouble( PyDict_GetItemString(self->params, "cyl_theta") );
+    self->model->background = PyFloat_AsDouble( PyDict_GetItemString(self->params, "background") );
+    self->model->radius = PyFloat_AsDouble( PyDict_GetItemString(self->params, "radius") );
+    self->model->contrast = PyFloat_AsDouble( PyDict_GetItemString(self->params, "contrast") );
+    self->model->cyl_phi = PyFloat_AsDouble( PyDict_GetItemString(self->params, "cyl_phi") );
+
+    // Read in dispersion parameters
+    DispersionVisitor* visitor = new DispersionVisitor();
+    PyObject* disp_dict = PyDict_New();
+    disp_dict = PyDict_GetItemString(self->dispersion, "radius");
+    self->model->radius.dispersion->accept_as_destination(visitor, self->model->radius.dispersion, disp_dict);
+
+    disp_dict = PyDict_GetItemString(self->dispersion, "length");
+    self->model->length.dispersion->accept_as_destination(visitor, self->model->length.dispersion, disp_dict);
+
+    disp_dict = PyDict_GetItemString(self->dispersion, "cyl_theta");
+    self->model->cyl_theta.dispersion->accept_as_destination(visitor, self->model->cyl_theta.dispersion, disp_dict);
+
+    disp_dict = PyDict_GetItemString(self->dispersion, "cyl_phi");
+    self->model->cyl_phi.dispersion->accept_as_destination(visitor, self->model->cyl_phi.dispersion, disp_dict);
+
+	// Get input and determine whether we have to supply a 1D or 2D return value.
+	if ( !PyArg_ParseTuple(args,"O",&pars) ) {
+	    PyErr_SetString(CCylinderModelError,
+	    	"CCylinderModel.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(CCylinderModelError,
+	    		"CCylinderModel.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 = CCylinderModel_readDouble(PyList_GET_ITEM(pars,0));
+	    qy_value = CCylinderModel_readDouble(PyList_GET_ITEM(pars,1));
+	    return Py_BuildValue("d",(*(self->model))(qx_value,qy_value));
+
+	} else {
+
+		// We have a scalar q, we will evaluate I(q)
+		qx_value = CCylinderModel_readDouble(pars);
+		return Py_BuildValue("d",(*(self->model))(qx_value));
+	}
+}
+
+static PyObject * reset(CCylinderModel *self, PyObject *args) {
+
+
+    return Py_BuildValue("d",0.0);
+}
+
+static PyObject * set_dispersion(CCylinderModel *self, PyObject *args) {
+	PyObject * disp;
+	const char * par_name;
+
+	if ( !PyArg_ParseTuple(args,"sO", &par_name, &disp) ) {
+	    PyErr_SetString(CCylinderModelError,
+	    	"CCylinderModel.set_dispersion expects a DispersionModel object.");
+		return NULL;
+	}
+	void *temp = PyCObject_AsVoidPtr(disp);
+	DispersionModel * dispersion = static_cast<DispersionModel *>(temp);
+
+
+	// Ugliness necessary to go from python to C
+	if (!strcmp(par_name, "radius")) {
+		self->model->radius.dispersion = dispersion;
+	} else if (!strcmp(par_name, "length")) {
+		self->model->length.dispersion = dispersion;
+	} else if (!strcmp(par_name, "cyl_theta")) {
+		self->model->cyl_theta.dispersion = dispersion;
+	} else if (!strcmp(par_name, "cyl_phi")) {
+		self->model->cyl_phi.dispersion = dispersion;
+	} else {
+	    PyErr_SetString(CCylinderModelError,
+	    	"CCylinderModel.set_dispersion expects a valid parameter name.");
+		return NULL;
+	}
+
+	DispersionVisitor* visitor = new DispersionVisitor();
+	PyObject * disp_dict = PyDict_New();
+	dispersion->accept_as_source(visitor, dispersion, disp_dict);
+	PyDict_SetItemString(self->dispersion, par_name, disp_dict);
+    return Py_BuildValue("i",1);
+}
+
+
+static PyMethodDef CCylinderModel_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"},
+    {"set_dispersion",      (PyCFunction)set_dispersion     , METH_VARARGS,
+      "Set the dispersion model for a given parameter"},
+   {NULL}
+};
+
+static PyTypeObject CCylinderModelType = {
+    PyObject_HEAD_INIT(NULL)
+    0,                         /*ob_size*/
+    "CCylinderModel",             /*tp_name*/
+    sizeof(CCylinderModel),             /*tp_basicsize*/
+    0,                         /*tp_itemsize*/
+    (destructor)CCylinderModel_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*/
+    "CCylinderModel objects",           /* tp_doc */
+    0,		               /* tp_traverse */
+    0,		               /* tp_clear */
+    0,		               /* tp_richcompare */
+    0,		               /* tp_weaklistoffset */
+    0,		               /* tp_iter */
+    0,		               /* tp_iternext */
+    CCylinderModel_methods,             /* tp_methods */
+    CCylinderModel_members,             /* tp_members */
+    0,                         /* tp_getset */
+    0,                         /* tp_base */
+    0,                         /* tp_dict */
+    0,                         /* tp_descr_get */
+    0,                         /* tp_descr_set */
+    0,                         /* tp_dictoffset */
+    (initproc)CCylinderModel_init,      /* tp_init */
+    0,                         /* tp_alloc */
+    CCylinderModel_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 addCCylinderModel(PyObject *module) {
+	PyObject *d;
+
+    if (PyType_Ready(&CCylinderModelType) < 0)
+        return;
+
+    Py_INCREF(&CCylinderModelType);
+    PyModule_AddObject(module, "CCylinderModel", (PyObject *)&CCylinderModelType);
+
+    d = PyModule_GetDict(module);
+    CCylinderModelError = PyErr_NewException("CCylinderModel.error", NULL, NULL);
+    PyDict_SetItemString(d, "CCylinderModelError", CCylinderModelError);
+}
+
Index: sansmodels/src/sans/models/c_models/CEllipsoidModel.cpp
===================================================================
--- sansmodels/src/sans/models/c_models/CEllipsoidModel.cpp	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
+++ sansmodels/src/sans/models/c_models/CEllipsoidModel.cpp	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
@@ -0,0 +1,399 @@
+/**
+	This software was developed by the University of Tennessee as part of the
+	Distributed Data Analysis of Neutron Scattering Experiments (DANSE)
+	project funded by the US National Science Foundation.
+
+	If you use DANSE applications to do scientific research that leads to
+	publication, we ask that you acknowledge the use of the software with the
+	following sentence:
+
+	"This work benefited from DANSE software developed under NSF award DMR-0520547."
+
+	copyright 2008, University of Tennessee
+ */
+
+/** CEllipsoidModel
+ *
+ * C extension 
+ *
+ * WARNING: THIS FILE WAS GENERATED BY WRAPPERGENERATOR.PY
+ *          DO NOT MODIFY THIS FILE, MODIFY ellipsoid.h
+ *          AND RE-RUN THE GENERATOR SCRIPT
+ *
+ */
+ 
+extern "C" {
+#include <Python.h>
+#include "structmember.h"
+#include <stdio.h>
+#include <stdlib.h>
+#include <math.h>
+#include <time.h>
+#include "ellipsoid.h"
+}
+
+#include "models.hh"
+#include "dispersion_visitor.hh"
+
+/// Error object for raised exceptions
+static PyObject * CEllipsoidModelError = NULL;
+
+
+// Class definition
+typedef struct {
+    PyObject_HEAD
+    /// Parameters
+    PyObject * params;
+    /// Dispersion parameters
+    PyObject * dispersion;
+    /// Underlying model object
+    EllipsoidModel * model;
+    /// Log for unit testing
+    PyObject * log;
+} CEllipsoidModel;
+
+
+static void
+CEllipsoidModel_dealloc(CEllipsoidModel* self)
+{
+    self->ob_type->tp_free((PyObject*)self);
+    
+
+}
+
+static PyObject *
+CEllipsoidModel_new(PyTypeObject *type, PyObject *args, PyObject *kwds)
+{
+    CEllipsoidModel *self;
+    
+    self = (CEllipsoidModel *)type->tp_alloc(type, 0);
+   
+    return (PyObject *)self;
+}
+
+static int
+CEllipsoidModel_init(CEllipsoidModel *self, PyObject *args, PyObject *kwds)
+{
+    if (self != NULL) {
+    	
+    	// Create parameters
+        self->params = PyDict_New();
+        self->dispersion = PyDict_New();
+        self->model = new EllipsoidModel();
+        
+        // Initialize parameter dictionary
+        PyDict_SetItemString(self->params,"scale",Py_BuildValue("d",1.000000));
+        PyDict_SetItemString(self->params,"axis_theta",Py_BuildValue("d",1.570000));
+        PyDict_SetItemString(self->params,"radius_b",Py_BuildValue("d",400.000000));
+        PyDict_SetItemString(self->params,"radius_a",Py_BuildValue("d",20.000000));
+        PyDict_SetItemString(self->params,"axis_phi",Py_BuildValue("d",0.000000));
+        PyDict_SetItemString(self->params,"background",Py_BuildValue("d",0.000000));
+        PyDict_SetItemString(self->params,"contrast",Py_BuildValue("d",0.000003));
+        // Initialize dispersion / averaging parameter dict
+        DispersionVisitor* visitor = new DispersionVisitor();
+        PyObject * disp_dict;
+        disp_dict = PyDict_New();
+        self->model->radius_a.dispersion->accept_as_source(visitor, self->model->radius_a.dispersion, disp_dict);
+        PyDict_SetItemString(self->dispersion, "radius_a", disp_dict);
+        disp_dict = PyDict_New();
+        self->model->radius_b.dispersion->accept_as_source(visitor, self->model->radius_b.dispersion, disp_dict);
+        PyDict_SetItemString(self->dispersion, "radius_b", disp_dict);
+        disp_dict = PyDict_New();
+        self->model->axis_theta.dispersion->accept_as_source(visitor, self->model->axis_theta.dispersion, disp_dict);
+        PyDict_SetItemString(self->dispersion, "axis_theta", disp_dict);
+        disp_dict = PyDict_New();
+        self->model->axis_phi.dispersion->accept_as_source(visitor, self->model->axis_phi.dispersion, disp_dict);
+        PyDict_SetItemString(self->dispersion, "axis_phi", disp_dict);
+
+
+         
+        // Create empty log
+        self->log = PyDict_New();
+        
+        
+
+    }
+    return 0;
+}
+
+static PyMemberDef CEllipsoidModel_members[] = {
+    {"params", T_OBJECT, offsetof(CEllipsoidModel, params), 0,
+     "Parameters"},
+	{"dispersion", T_OBJECT, offsetof(CEllipsoidModel, dispersion), 0,
+	  "Dispersion parameters"},     
+    {"log", T_OBJECT, offsetof(CEllipsoidModel, log), 0,
+     "Log"},
+    {NULL}  /* Sentinel */
+};
+
+/** Read double from PyObject
+    @param p PyObject
+    @return double
+*/
+double CEllipsoidModel_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(CEllipsoidModel *self, PyObject *args) {
+	double q_value, phi_value;
+	PyObject* pars;
+	int npars;
+	
+	// Get parameters
+	
+	    // Reader parameter dictionary
+    self->model->scale = PyFloat_AsDouble( PyDict_GetItemString(self->params, "scale") );
+    self->model->axis_theta = PyFloat_AsDouble( PyDict_GetItemString(self->params, "axis_theta") );
+    self->model->radius_b = PyFloat_AsDouble( PyDict_GetItemString(self->params, "radius_b") );
+    self->model->radius_a = PyFloat_AsDouble( PyDict_GetItemString(self->params, "radius_a") );
+    self->model->axis_phi = PyFloat_AsDouble( PyDict_GetItemString(self->params, "axis_phi") );
+    self->model->background = PyFloat_AsDouble( PyDict_GetItemString(self->params, "background") );
+    self->model->contrast = PyFloat_AsDouble( PyDict_GetItemString(self->params, "contrast") );
+    // Read in dispersion parameters
+    PyObject* disp_dict;
+    DispersionVisitor* visitor = new DispersionVisitor();
+    disp_dict = PyDict_GetItemString(self->dispersion, "radius_a");
+    self->model->radius_a.dispersion->accept_as_destination(visitor, self->model->radius_a.dispersion, disp_dict);
+    disp_dict = PyDict_GetItemString(self->dispersion, "radius_b");
+    self->model->radius_b.dispersion->accept_as_destination(visitor, self->model->radius_b.dispersion, disp_dict);
+    disp_dict = PyDict_GetItemString(self->dispersion, "axis_theta");
+    self->model->axis_theta.dispersion->accept_as_destination(visitor, self->model->axis_theta.dispersion, disp_dict);
+    disp_dict = PyDict_GetItemString(self->dispersion, "axis_phi");
+    self->model->axis_phi.dispersion->accept_as_destination(visitor, self->model->axis_phi.dispersion, disp_dict);
+
+	
+	// Get input and determine whether we have to supply a 1D or 2D return value.
+	if ( !PyArg_ParseTuple(args,"O",&pars) ) {
+	    PyErr_SetString(CEllipsoidModelError, 
+	    	"CEllipsoidModel.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(CEllipsoidModelError, 
+	    		"CEllipsoidModel.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 = CEllipsoidModel_readDouble(PyList_GET_ITEM(pars,0));
+	    phi_value = CEllipsoidModel_readDouble(PyList_GET_ITEM(pars,1));
+	    // Skip zero
+	    if (q_value==0) {
+	    	return Py_BuildValue("d",0.0);
+	    }
+		return Py_BuildValue("d",(*(self->model)).evaluate_rphi(q_value,phi_value));
+
+	} else {
+
+		// We have a scalar q, we will evaluate I(q)
+		q_value = CEllipsoidModel_readDouble(pars);		
+		
+		return Py_BuildValue("d",(*(self->model))(q_value));
+	}	
+}
+
+/**
+ * Function to call to evaluate model in cartesian coordinates
+ * @param args: input q or [qx, qy]]
+ * @return: function value
+ */
+static PyObject * runXY(CEllipsoidModel *self, PyObject *args) {
+	double qx_value, qy_value;
+	PyObject* pars;
+	int npars;
+	
+	// Get parameters
+	
+	    // Reader parameter dictionary
+    self->model->scale = PyFloat_AsDouble( PyDict_GetItemString(self->params, "scale") );
+    self->model->axis_theta = PyFloat_AsDouble( PyDict_GetItemString(self->params, "axis_theta") );
+    self->model->radius_b = PyFloat_AsDouble( PyDict_GetItemString(self->params, "radius_b") );
+    self->model->radius_a = PyFloat_AsDouble( PyDict_GetItemString(self->params, "radius_a") );
+    self->model->axis_phi = PyFloat_AsDouble( PyDict_GetItemString(self->params, "axis_phi") );
+    self->model->background = PyFloat_AsDouble( PyDict_GetItemString(self->params, "background") );
+    self->model->contrast = PyFloat_AsDouble( PyDict_GetItemString(self->params, "contrast") );
+    // Read in dispersion parameters
+    PyObject* disp_dict;
+    DispersionVisitor* visitor = new DispersionVisitor();
+    disp_dict = PyDict_GetItemString(self->dispersion, "radius_a");
+    self->model->radius_a.dispersion->accept_as_destination(visitor, self->model->radius_a.dispersion, disp_dict);
+    disp_dict = PyDict_GetItemString(self->dispersion, "radius_b");
+    self->model->radius_b.dispersion->accept_as_destination(visitor, self->model->radius_b.dispersion, disp_dict);
+    disp_dict = PyDict_GetItemString(self->dispersion, "axis_theta");
+    self->model->axis_theta.dispersion->accept_as_destination(visitor, self->model->axis_theta.dispersion, disp_dict);
+    disp_dict = PyDict_GetItemString(self->dispersion, "axis_phi");
+    self->model->axis_phi.dispersion->accept_as_destination(visitor, self->model->axis_phi.dispersion, disp_dict);
+
+	
+	// Get input and determine whether we have to supply a 1D or 2D return value.
+	if ( !PyArg_ParseTuple(args,"O",&pars) ) {
+	    PyErr_SetString(CEllipsoidModelError, 
+	    	"CEllipsoidModel.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(CEllipsoidModelError, 
+	    		"CEllipsoidModel.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 = CEllipsoidModel_readDouble(PyList_GET_ITEM(pars,0));
+	    qy_value = CEllipsoidModel_readDouble(PyList_GET_ITEM(pars,1));
+	    return Py_BuildValue("d",(*(self->model))(qx_value,qy_value));
+
+	} else {
+
+		// We have a scalar q, we will evaluate I(q)
+		qx_value = CEllipsoidModel_readDouble(pars);		
+		
+		return Py_BuildValue("d",(*(self->model))(qx_value));
+	}	
+}
+
+static PyObject * reset(CEllipsoidModel *self, PyObject *args) {
+    
+
+    return Py_BuildValue("d",0.0);
+}
+
+static PyObject * set_dispersion(CEllipsoidModel *self, PyObject *args) {
+	PyObject * disp;
+	const char * par_name;
+
+	if ( !PyArg_ParseTuple(args,"sO", &par_name, &disp) ) {
+	    PyErr_SetString(CEllipsoidModelError,
+	    	"CEllipsoidModel.set_dispersion expects a DispersionModel object.");
+		return NULL;
+	}
+	void *temp = PyCObject_AsVoidPtr(disp);
+	DispersionModel * dispersion = static_cast<DispersionModel *>(temp);
+
+
+	// Ugliness necessary to go from python to C
+	    // TODO: refactor this
+    if (!strcmp(par_name, "radius_a")) {
+        self->model->radius_a.dispersion = dispersion;
+    } else    if (!strcmp(par_name, "radius_b")) {
+        self->model->radius_b.dispersion = dispersion;
+    } else    if (!strcmp(par_name, "axis_theta")) {
+        self->model->axis_theta.dispersion = dispersion;
+    } else    if (!strcmp(par_name, "axis_phi")) {
+        self->model->axis_phi.dispersion = dispersion;
+    } else {
+	    PyErr_SetString(CEllipsoidModelError,
+	    	"CEllipsoidModel.set_dispersion expects a valid parameter name.");
+		return NULL;
+	}
+
+	DispersionVisitor* visitor = new DispersionVisitor();
+	PyObject * disp_dict = PyDict_New();
+	dispersion->accept_as_source(visitor, dispersion, disp_dict);
+	PyDict_SetItemString(self->dispersion, par_name, disp_dict);
+    return Py_BuildValue("i",1);
+}
+
+
+static PyMethodDef CEllipsoidModel_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"},
+    {"set_dispersion",      (PyCFunction)set_dispersion     , METH_VARARGS,
+      "Set the dispersion model for a given parameter"},
+   {NULL}
+};
+
+static PyTypeObject CEllipsoidModelType = {
+    PyObject_HEAD_INIT(NULL)
+    0,                         /*ob_size*/
+    "CEllipsoidModel",             /*tp_name*/
+    sizeof(CEllipsoidModel),             /*tp_basicsize*/
+    0,                         /*tp_itemsize*/
+    (destructor)CEllipsoidModel_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*/
+    "CEllipsoidModel objects",           /* tp_doc */
+    0,		               /* tp_traverse */
+    0,		               /* tp_clear */
+    0,		               /* tp_richcompare */
+    0,		               /* tp_weaklistoffset */
+    0,		               /* tp_iter */
+    0,		               /* tp_iternext */
+    CEllipsoidModel_methods,             /* tp_methods */
+    CEllipsoidModel_members,             /* tp_members */
+    0,                         /* tp_getset */
+    0,                         /* tp_base */
+    0,                         /* tp_dict */
+    0,                         /* tp_descr_get */
+    0,                         /* tp_descr_set */
+    0,                         /* tp_dictoffset */
+    (initproc)CEllipsoidModel_init,      /* tp_init */
+    0,                         /* tp_alloc */
+    CEllipsoidModel_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 addCEllipsoidModel(PyObject *module) {
+	PyObject *d;
+	
+    if (PyType_Ready(&CEllipsoidModelType) < 0)
+        return;
+
+    Py_INCREF(&CEllipsoidModelType);
+    PyModule_AddObject(module, "CEllipsoidModel", (PyObject *)&CEllipsoidModelType);
+    
+    d = PyModule_GetDict(module);
+    CEllipsoidModelError = PyErr_NewException("CEllipsoidModel.error", NULL, NULL);
+    PyDict_SetItemString(d, "CEllipsoidModelError", CEllipsoidModelError);
+}
+
Index: sansmodels/src/sans/models/c_models/CEllipticalCylinderModel.cpp
===================================================================
--- sansmodels/src/sans/models/c_models/CEllipticalCylinderModel.cpp	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
+++ sansmodels/src/sans/models/c_models/CEllipticalCylinderModel.cpp	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
@@ -0,0 +1,423 @@
+/**
+	This software was developed by the University of Tennessee as part of the
+	Distributed Data Analysis of Neutron Scattering Experiments (DANSE)
+	project funded by the US National Science Foundation.
+
+	If you use DANSE applications to do scientific research that leads to
+	publication, we ask that you acknowledge the use of the software with the
+	following sentence:
+
+	"This work benefited from DANSE software developed under NSF award DMR-0520547."
+
+	copyright 2008, University of Tennessee
+ */
+
+/** CEllipticalCylinderModel
+ *
+ * C extension 
+ *
+ * WARNING: THIS FILE WAS GENERATED BY WRAPPERGENERATOR.PY
+ *          DO NOT MODIFY THIS FILE, MODIFY elliptical_cylinder.h
+ *          AND RE-RUN THE GENERATOR SCRIPT
+ *
+ */
+ 
+extern "C" {
+#include <Python.h>
+#include "structmember.h"
+#include <stdio.h>
+#include <stdlib.h>
+#include <math.h>
+#include <time.h>
+#include "elliptical_cylinder.h"
+}
+
+#include "models.hh"
+#include "dispersion_visitor.hh"
+
+/// Error object for raised exceptions
+static PyObject * CEllipticalCylinderModelError = NULL;
+
+
+// Class definition
+typedef struct {
+    PyObject_HEAD
+    /// Parameters
+    PyObject * params;
+    /// Dispersion parameters
+    PyObject * dispersion;
+    /// Underlying model object
+    EllipticalCylinderModel * model;
+    /// Log for unit testing
+    PyObject * log;
+} CEllipticalCylinderModel;
+
+
+static void
+CEllipticalCylinderModel_dealloc(CEllipticalCylinderModel* self)
+{
+    self->ob_type->tp_free((PyObject*)self);
+    
+
+}
+
+static PyObject *
+CEllipticalCylinderModel_new(PyTypeObject *type, PyObject *args, PyObject *kwds)
+{
+    CEllipticalCylinderModel *self;
+    
+    self = (CEllipticalCylinderModel *)type->tp_alloc(type, 0);
+   
+    return (PyObject *)self;
+}
+
+static int
+CEllipticalCylinderModel_init(CEllipticalCylinderModel *self, PyObject *args, PyObject *kwds)
+{
+    if (self != NULL) {
+    	
+    	// Create parameters
+        self->params = PyDict_New();
+        self->dispersion = PyDict_New();
+        self->model = new EllipticalCylinderModel();
+        
+        // Initialize parameter dictionary
+        PyDict_SetItemString(self->params,"scale",Py_BuildValue("d",1.000000));
+        PyDict_SetItemString(self->params,"cyl_psi",Py_BuildValue("d",0.000000));
+        PyDict_SetItemString(self->params,"length",Py_BuildValue("d",400.000000));
+        PyDict_SetItemString(self->params,"r_minor",Py_BuildValue("d",20.000000));
+        PyDict_SetItemString(self->params,"cyl_theta",Py_BuildValue("d",1.570000));
+        PyDict_SetItemString(self->params,"background",Py_BuildValue("d",0.000000));
+        PyDict_SetItemString(self->params,"r_ratio",Py_BuildValue("d",1.500000));
+        PyDict_SetItemString(self->params,"contrast",Py_BuildValue("d",0.000003));
+        PyDict_SetItemString(self->params,"cyl_phi",Py_BuildValue("d",0.000000));
+        // Initialize dispersion / averaging parameter dict
+        DispersionVisitor* visitor = new DispersionVisitor();
+        PyObject * disp_dict;
+        disp_dict = PyDict_New();
+        self->model->r_minor.dispersion->accept_as_source(visitor, self->model->r_minor.dispersion, disp_dict);
+        PyDict_SetItemString(self->dispersion, "r_minor", disp_dict);
+        disp_dict = PyDict_New();
+        self->model->r_ratio.dispersion->accept_as_source(visitor, self->model->r_ratio.dispersion, disp_dict);
+        PyDict_SetItemString(self->dispersion, "r_ratio", disp_dict);
+        disp_dict = PyDict_New();
+        self->model->length.dispersion->accept_as_source(visitor, self->model->length.dispersion, disp_dict);
+        PyDict_SetItemString(self->dispersion, "length", disp_dict);
+        disp_dict = PyDict_New();
+        self->model->cyl_theta.dispersion->accept_as_source(visitor, self->model->cyl_theta.dispersion, disp_dict);
+        PyDict_SetItemString(self->dispersion, "cyl_theta", disp_dict);
+        disp_dict = PyDict_New();
+        self->model->cyl_phi.dispersion->accept_as_source(visitor, self->model->cyl_phi.dispersion, disp_dict);
+        PyDict_SetItemString(self->dispersion, "cyl_phi", disp_dict);
+        disp_dict = PyDict_New();
+        self->model->cyl_psi.dispersion->accept_as_source(visitor, self->model->cyl_psi.dispersion, disp_dict);
+        PyDict_SetItemString(self->dispersion, "cyl_psi", disp_dict);
+
+
+         
+        // Create empty log
+        self->log = PyDict_New();
+        
+        
+
+    }
+    return 0;
+}
+
+static PyMemberDef CEllipticalCylinderModel_members[] = {
+    {"params", T_OBJECT, offsetof(CEllipticalCylinderModel, params), 0,
+     "Parameters"},
+	{"dispersion", T_OBJECT, offsetof(CEllipticalCylinderModel, dispersion), 0,
+	  "Dispersion parameters"},     
+    {"log", T_OBJECT, offsetof(CEllipticalCylinderModel, log), 0,
+     "Log"},
+    {NULL}  /* Sentinel */
+};
+
+/** Read double from PyObject
+    @param p PyObject
+    @return double
+*/
+double CEllipticalCylinderModel_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(CEllipticalCylinderModel *self, PyObject *args) {
+	double q_value, phi_value;
+	PyObject* pars;
+	int npars;
+	
+	// Get parameters
+	
+	    // Reader parameter dictionary
+    self->model->scale = PyFloat_AsDouble( PyDict_GetItemString(self->params, "scale") );
+    self->model->cyl_psi = PyFloat_AsDouble( PyDict_GetItemString(self->params, "cyl_psi") );
+    self->model->length = PyFloat_AsDouble( PyDict_GetItemString(self->params, "length") );
+    self->model->r_minor = PyFloat_AsDouble( PyDict_GetItemString(self->params, "r_minor") );
+    self->model->cyl_theta = PyFloat_AsDouble( PyDict_GetItemString(self->params, "cyl_theta") );
+    self->model->background = PyFloat_AsDouble( PyDict_GetItemString(self->params, "background") );
+    self->model->r_ratio = PyFloat_AsDouble( PyDict_GetItemString(self->params, "r_ratio") );
+    self->model->contrast = PyFloat_AsDouble( PyDict_GetItemString(self->params, "contrast") );
+    self->model->cyl_phi = PyFloat_AsDouble( PyDict_GetItemString(self->params, "cyl_phi") );
+    // Read in dispersion parameters
+    PyObject* disp_dict;
+    DispersionVisitor* visitor = new DispersionVisitor();
+    disp_dict = PyDict_GetItemString(self->dispersion, "r_minor");
+    self->model->r_minor.dispersion->accept_as_destination(visitor, self->model->r_minor.dispersion, disp_dict);
+    disp_dict = PyDict_GetItemString(self->dispersion, "r_ratio");
+    self->model->r_ratio.dispersion->accept_as_destination(visitor, self->model->r_ratio.dispersion, disp_dict);
+    disp_dict = PyDict_GetItemString(self->dispersion, "length");
+    self->model->length.dispersion->accept_as_destination(visitor, self->model->length.dispersion, disp_dict);
+    disp_dict = PyDict_GetItemString(self->dispersion, "cyl_theta");
+    self->model->cyl_theta.dispersion->accept_as_destination(visitor, self->model->cyl_theta.dispersion, disp_dict);
+    disp_dict = PyDict_GetItemString(self->dispersion, "cyl_phi");
+    self->model->cyl_phi.dispersion->accept_as_destination(visitor, self->model->cyl_phi.dispersion, disp_dict);
+    disp_dict = PyDict_GetItemString(self->dispersion, "cyl_psi");
+    self->model->cyl_psi.dispersion->accept_as_destination(visitor, self->model->cyl_psi.dispersion, disp_dict);
+
+	
+	// Get input and determine whether we have to supply a 1D or 2D return value.
+	if ( !PyArg_ParseTuple(args,"O",&pars) ) {
+	    PyErr_SetString(CEllipticalCylinderModelError, 
+	    	"CEllipticalCylinderModel.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(CEllipticalCylinderModelError, 
+	    		"CEllipticalCylinderModel.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 = CEllipticalCylinderModel_readDouble(PyList_GET_ITEM(pars,0));
+	    phi_value = CEllipticalCylinderModel_readDouble(PyList_GET_ITEM(pars,1));
+	    // Skip zero
+	    if (q_value==0) {
+	    	return Py_BuildValue("d",0.0);
+	    }
+		return Py_BuildValue("d",(*(self->model)).evaluate_rphi(q_value,phi_value));
+
+	} else {
+
+		// We have a scalar q, we will evaluate I(q)
+		q_value = CEllipticalCylinderModel_readDouble(pars);		
+		
+		return Py_BuildValue("d",(*(self->model))(q_value));
+	}	
+}
+
+/**
+ * Function to call to evaluate model in cartesian coordinates
+ * @param args: input q or [qx, qy]]
+ * @return: function value
+ */
+static PyObject * runXY(CEllipticalCylinderModel *self, PyObject *args) {
+	double qx_value, qy_value;
+	PyObject* pars;
+	int npars;
+	
+	// Get parameters
+	
+	    // Reader parameter dictionary
+    self->model->scale = PyFloat_AsDouble( PyDict_GetItemString(self->params, "scale") );
+    self->model->cyl_psi = PyFloat_AsDouble( PyDict_GetItemString(self->params, "cyl_psi") );
+    self->model->length = PyFloat_AsDouble( PyDict_GetItemString(self->params, "length") );
+    self->model->r_minor = PyFloat_AsDouble( PyDict_GetItemString(self->params, "r_minor") );
+    self->model->cyl_theta = PyFloat_AsDouble( PyDict_GetItemString(self->params, "cyl_theta") );
+    self->model->background = PyFloat_AsDouble( PyDict_GetItemString(self->params, "background") );
+    self->model->r_ratio = PyFloat_AsDouble( PyDict_GetItemString(self->params, "r_ratio") );
+    self->model->contrast = PyFloat_AsDouble( PyDict_GetItemString(self->params, "contrast") );
+    self->model->cyl_phi = PyFloat_AsDouble( PyDict_GetItemString(self->params, "cyl_phi") );
+    // Read in dispersion parameters
+    PyObject* disp_dict;
+    DispersionVisitor* visitor = new DispersionVisitor();
+    disp_dict = PyDict_GetItemString(self->dispersion, "r_minor");
+    self->model->r_minor.dispersion->accept_as_destination(visitor, self->model->r_minor.dispersion, disp_dict);
+    disp_dict = PyDict_GetItemString(self->dispersion, "r_ratio");
+    self->model->r_ratio.dispersion->accept_as_destination(visitor, self->model->r_ratio.dispersion, disp_dict);
+    disp_dict = PyDict_GetItemString(self->dispersion, "length");
+    self->model->length.dispersion->accept_as_destination(visitor, self->model->length.dispersion, disp_dict);
+    disp_dict = PyDict_GetItemString(self->dispersion, "cyl_theta");
+    self->model->cyl_theta.dispersion->accept_as_destination(visitor, self->model->cyl_theta.dispersion, disp_dict);
+    disp_dict = PyDict_GetItemString(self->dispersion, "cyl_phi");
+    self->model->cyl_phi.dispersion->accept_as_destination(visitor, self->model->cyl_phi.dispersion, disp_dict);
+    disp_dict = PyDict_GetItemString(self->dispersion, "cyl_psi");
+    self->model->cyl_psi.dispersion->accept_as_destination(visitor, self->model->cyl_psi.dispersion, disp_dict);
+
+	
+	// Get input and determine whether we have to supply a 1D or 2D return value.
+	if ( !PyArg_ParseTuple(args,"O",&pars) ) {
+	    PyErr_SetString(CEllipticalCylinderModelError, 
+	    	"CEllipticalCylinderModel.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(CEllipticalCylinderModelError, 
+	    		"CEllipticalCylinderModel.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 = CEllipticalCylinderModel_readDouble(PyList_GET_ITEM(pars,0));
+	    qy_value = CEllipticalCylinderModel_readDouble(PyList_GET_ITEM(pars,1));
+	    return Py_BuildValue("d",(*(self->model))(qx_value,qy_value));
+
+	} else {
+
+		// We have a scalar q, we will evaluate I(q)
+		qx_value = CEllipticalCylinderModel_readDouble(pars);		
+		
+		return Py_BuildValue("d",(*(self->model))(qx_value));
+	}	
+}
+
+static PyObject * reset(CEllipticalCylinderModel *self, PyObject *args) {
+    
+
+    return Py_BuildValue("d",0.0);
+}
+
+static PyObject * set_dispersion(CEllipticalCylinderModel *self, PyObject *args) {
+	PyObject * disp;
+	const char * par_name;
+
+	if ( !PyArg_ParseTuple(args,"sO", &par_name, &disp) ) {
+	    PyErr_SetString(CEllipticalCylinderModelError,
+	    	"CEllipticalCylinderModel.set_dispersion expects a DispersionModel object.");
+		return NULL;
+	}
+	void *temp = PyCObject_AsVoidPtr(disp);
+	DispersionModel * dispersion = static_cast<DispersionModel *>(temp);
+
+
+	// Ugliness necessary to go from python to C
+	    // TODO: refactor this
+    if (!strcmp(par_name, "r_minor")) {
+        self->model->r_minor.dispersion = dispersion;
+    } else    if (!strcmp(par_name, "r_ratio")) {
+        self->model->r_ratio.dispersion = dispersion;
+    } else    if (!strcmp(par_name, "length")) {
+        self->model->length.dispersion = dispersion;
+    } else    if (!strcmp(par_name, "cyl_theta")) {
+        self->model->cyl_theta.dispersion = dispersion;
+    } else    if (!strcmp(par_name, "cyl_phi")) {
+        self->model->cyl_phi.dispersion = dispersion;
+    } else    if (!strcmp(par_name, "cyl_psi")) {
+        self->model->cyl_psi.dispersion = dispersion;
+    } else {
+	    PyErr_SetString(CEllipticalCylinderModelError,
+	    	"CEllipticalCylinderModel.set_dispersion expects a valid parameter name.");
+		return NULL;
+	}
+
+	DispersionVisitor* visitor = new DispersionVisitor();
+	PyObject * disp_dict = PyDict_New();
+	dispersion->accept_as_source(visitor, dispersion, disp_dict);
+	PyDict_SetItemString(self->dispersion, par_name, disp_dict);
+    return Py_BuildValue("i",1);
+}
+
+
+static PyMethodDef CEllipticalCylinderModel_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"},
+    {"set_dispersion",      (PyCFunction)set_dispersion     , METH_VARARGS,
+      "Set the dispersion model for a given parameter"},
+   {NULL}
+};
+
+static PyTypeObject CEllipticalCylinderModelType = {
+    PyObject_HEAD_INIT(NULL)
+    0,                         /*ob_size*/
+    "CEllipticalCylinderModel",             /*tp_name*/
+    sizeof(CEllipticalCylinderModel),             /*tp_basicsize*/
+    0,                         /*tp_itemsize*/
+    (destructor)CEllipticalCylinderModel_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*/
+    "CEllipticalCylinderModel objects",           /* tp_doc */
+    0,		               /* tp_traverse */
+    0,		               /* tp_clear */
+    0,		               /* tp_richcompare */
+    0,		               /* tp_weaklistoffset */
+    0,		               /* tp_iter */
+    0,		               /* tp_iternext */
+    CEllipticalCylinderModel_methods,             /* tp_methods */
+    CEllipticalCylinderModel_members,             /* tp_members */
+    0,                         /* tp_getset */
+    0,                         /* tp_base */
+    0,                         /* tp_dict */
+    0,                         /* tp_descr_get */
+    0,                         /* tp_descr_set */
+    0,                         /* tp_dictoffset */
+    (initproc)CEllipticalCylinderModel_init,      /* tp_init */
+    0,                         /* tp_alloc */
+    CEllipticalCylinderModel_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 addCEllipticalCylinderModel(PyObject *module) {
+	PyObject *d;
+	
+    if (PyType_Ready(&CEllipticalCylinderModelType) < 0)
+        return;
+
+    Py_INCREF(&CEllipticalCylinderModelType);
+    PyModule_AddObject(module, "CEllipticalCylinderModel", (PyObject *)&CEllipticalCylinderModelType);
+    
+    d = PyModule_GetDict(module);
+    CEllipticalCylinderModelError = PyErr_NewException("CEllipticalCylinderModel.error", NULL, NULL);
+    PyDict_SetItemString(d, "CEllipticalCylinderModelError", CEllipticalCylinderModelError);
+}
+
Index: sansmodels/src/sans/models/c_models/CSphereModel.cpp
===================================================================
--- sansmodels/src/sans/models/c_models/CSphereModel.cpp	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
+++ sansmodels/src/sans/models/c_models/CSphereModel.cpp	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
@@ -0,0 +1,363 @@
+/**
+	This software was developed by the University of Tennessee as part of the
+	Distributed Data Analysis of Neutron Scattering Experiments (DANSE)
+	project funded by the US National Science Foundation.
+
+	If you use DANSE applications to do scientific research that leads to
+	publication, we ask that you acknowledge the use of the software with the
+	following sentence:
+
+	"This work benefited from DANSE software developed under NSF award DMR-0520547."
+
+	copyright 2008, University of Tennessee
+ */
+
+/** CSphereModel
+ *
+ * C extension 
+ *
+ * WARNING: THIS FILE WAS GENERATED BY WRAPPERGENERATOR.PY
+ *          DO NOT MODIFY THIS FILE, MODIFY sphere.h
+ *          AND RE-RUN THE GENERATOR SCRIPT
+ *
+ */
+ 
+extern "C" {
+#include <Python.h>
+#include "structmember.h"
+#include <stdio.h>
+#include <stdlib.h>
+#include <math.h>
+#include <time.h>
+#include "sphere.h"
+}
+
+#include "models.hh"
+#include "dispersion_visitor.hh"
+
+/// Error object for raised exceptions
+static PyObject * CSphereModelError = NULL;
+
+
+// Class definition
+typedef struct {
+    PyObject_HEAD
+    /// Parameters
+    PyObject * params;
+    /// Dispersion parameters
+    PyObject * dispersion;
+    /// Underlying model object
+    SphereModel * model;
+    /// Log for unit testing
+    PyObject * log;
+} CSphereModel;
+
+
+static void
+CSphereModel_dealloc(CSphereModel* self)
+{
+    self->ob_type->tp_free((PyObject*)self);
+    
+
+}
+
+static PyObject *
+CSphereModel_new(PyTypeObject *type, PyObject *args, PyObject *kwds)
+{
+    CSphereModel *self;
+    
+    self = (CSphereModel *)type->tp_alloc(type, 0);
+   
+    return (PyObject *)self;
+}
+
+static int
+CSphereModel_init(CSphereModel *self, PyObject *args, PyObject *kwds)
+{
+    if (self != NULL) {
+    	
+    	// Create parameters
+        self->params = PyDict_New();
+        self->dispersion = PyDict_New();
+        self->model = new SphereModel();
+        
+        // Initialize parameter dictionary
+        PyDict_SetItemString(self->params,"scale",Py_BuildValue("d",0.000001));
+        PyDict_SetItemString(self->params,"radius",Py_BuildValue("d",60.000000));
+        PyDict_SetItemString(self->params,"background",Py_BuildValue("d",0.000000));
+        PyDict_SetItemString(self->params,"contrast",Py_BuildValue("d",1.000000));
+        // Initialize dispersion / averaging parameter dict
+        DispersionVisitor* visitor = new DispersionVisitor();
+        PyObject * disp_dict;
+        disp_dict = PyDict_New();
+        self->model->radius.dispersion->accept_as_source(visitor, self->model->radius.dispersion, disp_dict);
+        PyDict_SetItemString(self->dispersion, "radius", disp_dict);
+
+
+         
+        // Create empty log
+        self->log = PyDict_New();
+        
+        
+
+    }
+    return 0;
+}
+
+static PyMemberDef CSphereModel_members[] = {
+    {"params", T_OBJECT, offsetof(CSphereModel, params), 0,
+     "Parameters"},
+	{"dispersion", T_OBJECT, offsetof(CSphereModel, dispersion), 0,
+	  "Dispersion parameters"},     
+    {"log", T_OBJECT, offsetof(CSphereModel, log), 0,
+     "Log"},
+    {NULL}  /* Sentinel */
+};
+
+/** Read double from PyObject
+    @param p PyObject
+    @return double
+*/
+double CSphereModel_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(CSphereModel *self, PyObject *args) {
+	double q_value, phi_value;
+	PyObject* pars;
+	int npars;
+	
+	// Get parameters
+	
+	    // Reader parameter dictionary
+    self->model->scale = PyFloat_AsDouble( PyDict_GetItemString(self->params, "scale") );
+    self->model->radius = PyFloat_AsDouble( PyDict_GetItemString(self->params, "radius") );
+    self->model->background = PyFloat_AsDouble( PyDict_GetItemString(self->params, "background") );
+    self->model->contrast = PyFloat_AsDouble( PyDict_GetItemString(self->params, "contrast") );
+    // Read in dispersion parameters
+    PyObject* disp_dict;
+    DispersionVisitor* visitor = new DispersionVisitor();
+    disp_dict = PyDict_GetItemString(self->dispersion, "radius");
+    self->model->radius.dispersion->accept_as_destination(visitor, self->model->radius.dispersion, disp_dict);
+
+	
+	// Get input and determine whether we have to supply a 1D or 2D return value.
+	if ( !PyArg_ParseTuple(args,"O",&pars) ) {
+	    PyErr_SetString(CSphereModelError, 
+	    	"CSphereModel.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(CSphereModelError, 
+	    		"CSphereModel.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 = CSphereModel_readDouble(PyList_GET_ITEM(pars,0));
+	    phi_value = CSphereModel_readDouble(PyList_GET_ITEM(pars,1));
+	    // Skip zero
+	    if (q_value==0) {
+	    	return Py_BuildValue("d",0.0);
+	    }
+		return Py_BuildValue("d",(*(self->model)).evaluate_rphi(q_value,phi_value));
+
+	} else {
+
+		// We have a scalar q, we will evaluate I(q)
+		q_value = CSphereModel_readDouble(pars);		
+		
+		return Py_BuildValue("d",(*(self->model))(q_value));
+	}	
+}
+
+/**
+ * Function to call to evaluate model in cartesian coordinates
+ * @param args: input q or [qx, qy]]
+ * @return: function value
+ */
+static PyObject * runXY(CSphereModel *self, PyObject *args) {
+	double qx_value, qy_value;
+	PyObject* pars;
+	int npars;
+	
+	// Get parameters
+	
+	    // Reader parameter dictionary
+    self->model->scale = PyFloat_AsDouble( PyDict_GetItemString(self->params, "scale") );
+    self->model->radius = PyFloat_AsDouble( PyDict_GetItemString(self->params, "radius") );
+    self->model->background = PyFloat_AsDouble( PyDict_GetItemString(self->params, "background") );
+    self->model->contrast = PyFloat_AsDouble( PyDict_GetItemString(self->params, "contrast") );
+    // Read in dispersion parameters
+    PyObject* disp_dict;
+    DispersionVisitor* visitor = new DispersionVisitor();
+    disp_dict = PyDict_GetItemString(self->dispersion, "radius");
+    self->model->radius.dispersion->accept_as_destination(visitor, self->model->radius.dispersion, disp_dict);
+
+	
+	// Get input and determine whether we have to supply a 1D or 2D return value.
+	if ( !PyArg_ParseTuple(args,"O",&pars) ) {
+	    PyErr_SetString(CSphereModelError, 
+	    	"CSphereModel.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(CSphereModelError, 
+	    		"CSphereModel.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 = CSphereModel_readDouble(PyList_GET_ITEM(pars,0));
+	    qy_value = CSphereModel_readDouble(PyList_GET_ITEM(pars,1));
+	    return Py_BuildValue("d",(*(self->model))(qx_value,qy_value));
+
+	} else {
+
+		// We have a scalar q, we will evaluate I(q)
+		qx_value = CSphereModel_readDouble(pars);		
+		
+		return Py_BuildValue("d",(*(self->model))(qx_value));
+	}	
+}
+
+static PyObject * reset(CSphereModel *self, PyObject *args) {
+    
+
+    return Py_BuildValue("d",0.0);
+}
+
+static PyObject * set_dispersion(CSphereModel *self, PyObject *args) {
+	PyObject * disp;
+	const char * par_name;
+
+	if ( !PyArg_ParseTuple(args,"sO", &par_name, &disp) ) {
+	    PyErr_SetString(CSphereModelError,
+	    	"CSphereModel.set_dispersion expects a DispersionModel object.");
+		return NULL;
+	}
+	void *temp = PyCObject_AsVoidPtr(disp);
+	DispersionModel * dispersion = static_cast<DispersionModel *>(temp);
+
+
+	// Ugliness necessary to go from python to C
+	    // TODO: refactor this
+    if (!strcmp(par_name, "radius")) {
+        self->model->radius.dispersion = dispersion;
+    } else {
+	    PyErr_SetString(CSphereModelError,
+	    	"CSphereModel.set_dispersion expects a valid parameter name.");
+		return NULL;
+	}
+
+	DispersionVisitor* visitor = new DispersionVisitor();
+	PyObject * disp_dict = PyDict_New();
+	dispersion->accept_as_source(visitor, dispersion, disp_dict);
+	PyDict_SetItemString(self->dispersion, par_name, disp_dict);
+    return Py_BuildValue("i",1);
+}
+
+
+static PyMethodDef CSphereModel_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"},
+    {"set_dispersion",      (PyCFunction)set_dispersion     , METH_VARARGS,
+      "Set the dispersion model for a given parameter"},
+   {NULL}
+};
+
+static PyTypeObject CSphereModelType = {
+    PyObject_HEAD_INIT(NULL)
+    0,                         /*ob_size*/
+    "CSphereModel",             /*tp_name*/
+    sizeof(CSphereModel),             /*tp_basicsize*/
+    0,                         /*tp_itemsize*/
+    (destructor)CSphereModel_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*/
+    "CSphereModel objects",           /* tp_doc */
+    0,		               /* tp_traverse */
+    0,		               /* tp_clear */
+    0,		               /* tp_richcompare */
+    0,		               /* tp_weaklistoffset */
+    0,		               /* tp_iter */
+    0,		               /* tp_iternext */
+    CSphereModel_methods,             /* tp_methods */
+    CSphereModel_members,             /* tp_members */
+    0,                         /* tp_getset */
+    0,                         /* tp_base */
+    0,                         /* tp_dict */
+    0,                         /* tp_descr_get */
+    0,                         /* tp_descr_set */
+    0,                         /* tp_dictoffset */
+    (initproc)CSphereModel_init,      /* tp_init */
+    0,                         /* tp_alloc */
+    CSphereModel_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 addCSphereModel(PyObject *module) {
+	PyObject *d;
+	
+    if (PyType_Ready(&CSphereModelType) < 0)
+        return;
+
+    Py_INCREF(&CSphereModelType);
+    PyModule_AddObject(module, "CSphereModel", (PyObject *)&CSphereModelType);
+    
+    d = PyModule_GetDict(module);
+    CSphereModelError = PyErr_NewException("CSphereModel.error", NULL, NULL);
+    PyDict_SetItemString(d, "CSphereModelError", CSphereModelError);
+}
+
Index: sansmodels/src/sans/models/c_models/WrapperGenerator.py
===================================================================
--- sansmodels/src/sans/models/c_models/WrapperGenerator.py	(revision af03ddd0932915d6659404ab7824e4223ea2944f)
+++ sansmodels/src/sans/models/c_models/WrapperGenerator.py	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
@@ -213,6 +213,6 @@
             param_str = "// Initialize parameter dictionary\n"            
             for par in self.params:
-                param_str += "        PyDict_SetItemString(self->params,\"%s\",Py_BuildValue(\"d\", self->model->%s()));\n" % \
-                    (par, par)
+                param_str += "        PyDict_SetItemString(self->params,\"%s\",Py_BuildValue(\"d\",%f));\n" % \
+                    (par, self.params[par])
 
             param_str += "        // Initialize dispersion / averaging parameter dict\n"
Index: sansmodels/src/sans/models/c_models/coreshellcylinder.cpp
===================================================================
--- sansmodels/src/sans/models/c_models/coreshellcylinder.cpp	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
+++ sansmodels/src/sans/models/c_models/coreshellcylinder.cpp	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
@@ -0,0 +1,214 @@
+/**
+	This software was developed by the University of Tennessee as part of the
+	Distributed Data Analysis of Neutron Scattering Experiments (DANSE)
+	project funded by the US National Science Foundation.
+
+	If you use DANSE applications to do scientific research that leads to
+	publication, we ask that you acknowledge the use of the software with the
+	following sentence:
+
+	"This work benefited from DANSE software developed under NSF award DMR-0520547."
+
+	copyright 2008, University of Tennessee
+ */
+
+/**
+ * Scattering model classes
+ * The classes use the IGOR library found in
+ *   sansmodels/src/libigor
+ *
+ *	TODO: refactor so that we pull in the old sansmodels.c_extensions
+ */
+
+#include <math.h>
+#include "models.hh"
+#include "parameters.hh"
+#include <stdio.h>
+using namespace std;
+
+extern "C" {
+	#include "libCylinder.h"
+	#include "core_shell_cylinder.h"
+}
+
+CoreShellCylinderModel :: CoreShellCylinderModel() {
+	scale      = Parameter(1.0);
+	radius     = Parameter(20.0, true);
+	radius.set_min(0.0);
+	thickness  = Parameter(10.0, true);
+	thickness.set_min(0.0);
+	length     = Parameter(400.0, true);
+	length.set_min(0.0);
+	core_sld   = Parameter(1.e-6);
+	shell_sld  = Parameter(4.e-6);
+	solvent_sld= Parameter(1.e-6);
+	background = Parameter(0.0);
+	axis_theta = Parameter(0.0, true);
+	axis_phi   = Parameter(0.0, true);
+}
+
+/**
+ * Function to evaluate 1D scattering function
+ * The NIST IGOR library is used for the actual calculation.
+ * @param q: q-value
+ * @return: function value
+ */
+double CoreShellCylinderModel :: operator()(double q) {
+	double dp[8];
+
+	dp[0] = scale();
+	dp[1] = radius();
+	dp[2] = thickness();
+	dp[3] = length();
+	dp[4] = core_sld();
+	dp[5] = shell_sld();
+	dp[6] = solvent_sld();
+	dp[7] = 0.0;
+
+	// Get the dispersion points for the radius
+	vector<WeightPoint> weights_rad;
+	radius.get_weights(weights_rad);
+
+	// Get the dispersion points for the thickness
+	vector<WeightPoint> weights_thick;
+	thickness.get_weights(weights_thick);
+
+	// Get the dispersion points for the length
+	vector<WeightPoint> weights_len;
+	length.get_weights(weights_len);
+
+	// Perform the computation, with all weight points
+	double sum = 0.0;
+	double norm = 0.0;
+
+	// Loop over radius weight points
+	for(int i=0; i<weights_rad.size(); i++) {
+		dp[1] = weights_rad[i].value;
+
+		// Loop over length weight points
+		for(int j=0; j<weights_len.size(); j++) {
+			dp[3] = weights_len[j].value;
+
+			// Loop over thickness weight points
+			for(int k=0; k<weights_thick.size(); k++) {
+				dp[2] = weights_thick[k].value;
+
+				sum += weights_rad[i].weight
+					* weights_len[j].weight
+					* weights_thick[k].weight
+					* CoreShellCylinder(dp, q);
+				norm += weights_rad[i].weight
+				* weights_len[j].weight
+				* weights_thick[k].weight;
+			}
+		}
+	}
+	return sum/norm + background();
+}
+
+/**
+ * Function to evaluate 2D scattering function
+ * @param q_x: value of Q along x
+ * @param q_y: value of Q along y
+ * @return: function value
+ */
+double CoreShellCylinderModel :: operator()(double qx, double qy) {
+	CoreShellCylinderParameters dp;
+	// Fill parameter array
+	dp.scale      = scale();
+	dp.radius     = radius();
+	dp.thickness  = thickness();
+	dp.length     = length();
+	dp.core_sld   = core_sld();
+	dp.shell_sld  = shell_sld();
+	dp.solvent_sld= solvent_sld();
+	dp.background = 0.0;
+	dp.axis_theta = axis_theta();
+	dp.axis_phi   = axis_phi();
+
+	// Get the dispersion points for the radius
+	vector<WeightPoint> weights_rad;
+	radius.get_weights(weights_rad);
+
+	// Get the dispersion points for the thickness
+	vector<WeightPoint> weights_thick;
+	thickness.get_weights(weights_thick);
+
+	// Get the dispersion points for the length
+	vector<WeightPoint> weights_len;
+	length.get_weights(weights_len);
+
+	// Get angular averaging for theta
+	vector<WeightPoint> weights_theta;
+	axis_theta.get_weights(weights_theta);
+
+	// Get angular averaging for phi
+	vector<WeightPoint> weights_phi;
+	axis_phi.get_weights(weights_phi);
+
+	// Perform the computation, with all weight points
+	double sum = 0.0;
+	double norm = 0.0;
+
+	// Loop over radius weight points
+	for(int i=0; i<weights_rad.size(); i++) {
+		dp.radius = weights_rad[i].value;
+
+
+		// Loop over length weight points
+		for(int j=0; j<weights_len.size(); j++) {
+			dp.length = weights_len[j].value;
+
+			// Loop over thickness weight points
+			for(int m=0; m<weights_thick.size(); m++) {
+				dp.thickness = weights_thick[m].value;
+
+			// Average over theta distribution
+			for(int k=0; k<weights_theta.size(); k++) {
+				dp.axis_theta = weights_theta[k].value;
+
+				// Average over phi distribution
+				for(int l=0; l<weights_phi.size(); l++) {
+					dp.axis_phi = weights_phi[l].value;
+
+					double _ptvalue = weights_rad[i].weight
+						* weights_len[j].weight
+						* weights_thick[m].weight
+						* weights_theta[k].weight
+						* weights_phi[l].weight
+						* core_shell_cylinder_analytical_2DXY(&dp, qx, qy);
+					if (weights_theta.size()>1) {
+						_ptvalue *= sin(weights_theta[k].value);
+					}
+					sum += _ptvalue;
+
+					norm += weights_rad[i].weight
+						* weights_len[j].weight
+						* weights_thick[m].weight
+						* weights_theta[k].weight
+						* weights_phi[l].weight;
+
+				}
+			}
+			}
+		}
+	}
+	// Averaging in theta needs an extra normalization
+	// factor to account for the sin(theta) term in the
+	// integration (see documentation).
+	if (weights_theta.size()>1) norm = norm / asin(1.0);
+	return sum/norm + background();
+}
+
+/**
+ * Function to evaluate 2D scattering function
+ * @param pars: parameters of the cylinder
+ * @param q: q-value
+ * @param phi: angle phi
+ * @return: function value
+ */
+double CoreShellCylinderModel :: evaluate_rphi(double q, double phi) {
+	double qx = q*cos(phi);
+	double qy = q*sin(phi);
+	return (*this).operator()(qx, qy);
+}
Index: sansmodels/src/sans/models/c_models/coreshellsphere.cpp
===================================================================
--- sansmodels/src/sans/models/c_models/coreshellsphere.cpp	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
+++ sansmodels/src/sans/models/c_models/coreshellsphere.cpp	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
@@ -0,0 +1,114 @@
+/**
+	This software was developed by the University of Tennessee as part of the
+	Distributed Data Analysis of Neutron Scattering Experiments (DANSE)
+	project funded by the US National Science Foundation.
+
+	If you use DANSE applications to do scientific research that leads to
+	publication, we ask that you acknowledge the use of the software with the
+	following sentence:
+
+	"This work benefited from DANSE software developed under NSF award DMR-0520547."
+
+	copyright 2008, University of Tennessee
+ */
+
+/**
+ * Scattering model classes
+ * The classes use the IGOR library found in
+ *   sansmodels/src/libigor
+ *
+ */
+
+#include <math.h>
+#include "models.hh"
+#include "parameters.hh"
+#include <stdio.h>
+using namespace std;
+
+extern "C" {
+	#include "libSphere.h"
+}
+
+CoreShellModel :: CoreShellModel() {
+	scale      = Parameter(1.0);
+	radius     = Parameter(60.0, true);
+	radius.set_min(0.0);
+	thickness  = Parameter(10.0, true);
+	thickness.set_min(0.0);
+	core_sld   = Parameter(1.e-6);
+	shell_sld  = Parameter(2.e-6);
+	solvent_sld = Parameter(3.e-6);
+	background = Parameter(0.0);
+}
+
+/**
+ * Function to evaluate 1D scattering function
+ * The NIST IGOR library is used for the actual calculation.
+ * @param q: q-value
+ * @return: function value
+ */
+double CoreShellModel :: operator()(double q) {
+	double dp[7];
+
+	// Fill parameter array for IGOR library
+	// Add the background after averaging
+
+	dp[0] = scale();
+	dp[1] = radius();
+	dp[2] = thickness();
+	dp[3] = core_sld();
+	dp[4] = shell_sld();
+	dp[5] = solvent_sld();
+	dp[6] = 0.0;
+
+
+	// Get the dispersion points for the radius
+	vector<WeightPoint> weights_rad;
+	radius.get_weights(weights_rad);
+
+	// Get the dispersion points for the thickness
+	vector<WeightPoint> weights_thick;
+	thickness.get_weights(weights_thick);
+
+	// Perform the computation, with all weight points
+	double sum = 0.0;
+	double norm = 0.0;
+
+	// Loop over radius weight points
+	for(int i=0; i<weights_rad.size(); i++) {
+		dp[1] = weights_rad[i].value;
+
+		// Loop over thickness weight points
+		for(int j=0; j<weights_thick.size(); j++) {
+			dp[2] = weights_thick[j].value;
+
+			sum += weights_rad[i].weight
+				* weights_thick[j].weight * CoreShellForm(dp, q);
+			norm += weights_rad[i].weight
+				* weights_thick[j].weight;
+		}
+	}
+	return sum/norm + background();
+}
+
+/**
+ * Function to evaluate 2D scattering function
+ * @param q_x: value of Q along x
+ * @param q_y: value of Q along y
+ * @return: function value
+ */
+double CoreShellModel :: operator()(double qx, double qy) {
+	double q = sqrt(qx*qx + qy*qy);
+	return (*this).operator()(q);
+}
+
+/**
+ * Function to evaluate 2D scattering function
+ * @param pars: parameters of the sphere
+ * @param q: q-value
+ * @param phi: angle phi
+ * @return: function value
+ */
+double CoreShellModel :: evaluate_rphi(double q, double phi) {
+	return (*this).operator()(q);
+}
Index: sansmodels/src/sans/models/c_models/cylinder.cpp
===================================================================
--- sansmodels/src/sans/models/c_models/cylinder.cpp	(revision af03ddd0932915d6659404ab7824e4223ea2944f)
+++ sansmodels/src/sans/models/c_models/cylinder.cpp	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
@@ -189,4 +189,5 @@
 
 	printf("Length = %g\n", c.length());
+	printf("I(Qx=%g,Qy=%g) = %g\n", 0.001, 0.001, c(0.001, 0.001));
 	printf("I(Q=%g) = %g\n", 0.001, c(0.001));
 	c.radius.dispersion = new GaussianDispersion();
Index: sansmodels/src/sans/models/c_models/dispersion_models.py
===================================================================
--- sansmodels/src/sans/models/c_models/dispersion_models.py	(revision af03ddd0932915d6659404ab7824e4223ea2944f)
+++ 	(revision )
@@ -1,100 +1,0 @@
-"""
-    This software was developed by the University of Tennessee as part of the
-    Distributed Data Analysis of Neutron Scattering Experiments (DANSE)
-    project funded by the US National Science Foundation.
-
-    If you use DANSE applications to do scientific research that leads to
-    publication, we ask that you acknowledge the use of the software with the
-    following sentence:
-
-    "This work benefited from DANSE software developed under NSF award DMR-0520547."
-
-    copyright 2008, University of Tennessee
-"""
-
-"""
-    Class definitions for python dispersion model for 
-    model parameters. These classes are bridges to the C++
-    dispersion object.
-    
-    The ArrayDispersion class takes in numpy arrays only.
-    
-    Usage:
-    These classes can be used to set the dispersion model of a SANS model
-    parameter:
-    
-        cyl = CylinderModel()
-        cyl.set_dispersion('radius', GaussianDispersion())
-    
-    
-    After the dispersion model is set, you can access it's
-    parameter through the dispersion dictionary:
-    
-        cyl.dispersion['radius']['width'] = 5.0
-    
-    TODO: For backward compatibility, the model parameters are still kept in
-    a dictionary. The next iteration of refactoring work should involve moving
-    away from value-based parameters to object-based parameter. We want to
-    store parameters as objects so that we can unify the 'params' and 'dispersion'
-    dictionaries into a single dictionary of parameter objects that hold the
-    complete information about the parameter (units, limits, dispersion model, etc...).  
-    
-    
-"""
-import sans_extension.c_models as c_models 
-
-class DispersionModel:
-    """
-        Python bridge class for a basic dispersion model 
-        class with a constant parameter value distribution
-    """
-    def __init__(self):
-        self.cdisp = c_models.new_dispersion_model()
-        
-    def set_weights(self, values, weights):
-        """
-            Set the weights of an array dispersion
-        """
-        message = "set_weights is not available for DispersionModel.\n"
-        message += "  Solution: Use an ArrayDispersion object"
-        raise "RuntimeError", message
-        
-class GaussianDispersion(DispersionModel):
-    """
-        Python bridge class for a dispersion model based 
-        on a Gaussian distribution.
-    """
-    def __init__(self):
-        self.cdisp = c_models.new_gaussian_model()
-        
-    def set_weights(self, values, weights):
-        """
-            Set the weights of an array dispersion
-        """
-        message = "set_weights is not available for GaussiantDispersion.\n"
-        message += "  Solution: Use an ArrayDispersion object"
-        raise "RuntimeError", message
-        
-class ArrayDispersion(DispersionModel):
-    """
-        Python bridge class for a dispersion model based on arrays.
-        The user has to set a weight distribution that
-        will be used in the averaging the model parameter
-        it is applied to. 
-    """
-    def __init__(self):
-        self.cdisp = c_models.new_array_model()
-        
-    def set_weights(self, values, weights):
-        """
-            Set the weights of an array dispersion
-            Only accept numpy arrays.
-            @param values: numpy array of values
-            @param weights: numpy array of weights for each value entry
-        """
-        if len(values) != len(weights):
-            raise ValueError, "ArrayDispersion.set_weights: given arrays are of different lengths"
-        
-        c_models.set_dispersion_weights(self.cdisp, values, weights)
-        
-        
Index: sansmodels/src/sans/models/c_models/ellipsoid.cpp
===================================================================
--- sansmodels/src/sans/models/c_models/ellipsoid.cpp	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
+++ sansmodels/src/sans/models/c_models/ellipsoid.cpp	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
@@ -0,0 +1,183 @@
+/**
+	This software was developed by the University of Tennessee as part of the
+	Distributed Data Analysis of Neutron Scattering Experiments (DANSE)
+	project funded by the US National Science Foundation.
+
+	If you use DANSE applications to do scientific research that leads to
+	publication, we ask that you acknowledge the use of the software with the
+	following sentence:
+
+	"This work benefited from DANSE software developed under NSF award DMR-0520547."
+
+	copyright 2008, University of Tennessee
+ */
+
+/**
+ * Scattering model classes
+ * The classes use the IGOR library found in
+ *   sansmodels/src/libigor
+ *
+ *	TODO: refactor so that we pull in the old sansmodels.c_extensions
+ */
+
+#include <math.h>
+#include "models.hh"
+#include "parameters.hh"
+#include <stdio.h>
+using namespace std;
+
+extern "C" {
+	#include "libCylinder.h"
+	#include "ellipsoid.h"
+}
+
+EllipsoidModel :: EllipsoidModel() {
+	scale      = Parameter(1.0);
+	radius_a   = Parameter(20.0, true);
+	radius_a.set_min(0.0);
+	radius_b   = Parameter(400.0, true);
+	radius_b.set_min(0.0);
+	contrast   = Parameter(3.e-6);
+	background = Parameter(0.0);
+	axis_theta  = Parameter(1.57, true);
+	axis_phi    = Parameter(0.0, true);
+}
+
+/**
+ * Function to evaluate 1D scattering function
+ * The NIST IGOR library is used for the actual calculation.
+ * @param q: q-value
+ * @return: function value
+ */
+double EllipsoidModel :: operator()(double q) {
+	double dp[5];
+
+	// Fill parameter array for IGOR library
+	// Add the background after averaging
+	dp[0] = scale();
+	dp[1] = radius_a();
+	dp[2] = radius_b();
+	dp[3] = contrast();
+	dp[4] = 0.0;
+
+	// Get the dispersion points for the radius_a
+	vector<WeightPoint> weights_rad_a;
+	radius_a.get_weights(weights_rad_a);
+
+	// Get the dispersion points for the radius_b
+	vector<WeightPoint> weights_rad_b;
+	radius_b.get_weights(weights_rad_b);
+
+	// Perform the computation, with all weight points
+	double sum = 0.0;
+	double norm = 0.0;
+
+	// Loop over radius_a weight points
+	for(int i=0; i<weights_rad_a.size(); i++) {
+		dp[1] = weights_rad_a[i].value;
+
+		// Loop over radius_b weight points
+		for(int j=0; j<weights_rad_b.size(); j++) {
+			dp[2] = weights_rad_b[j].value;
+
+			sum += weights_rad_a[i].weight
+				* weights_rad_b[j].weight * EllipsoidForm(dp, q);
+			norm += weights_rad_a[i].weight
+				* weights_rad_b[j].weight;
+		}
+	}
+	return sum/norm + background();
+}
+
+/**
+ * Function to evaluate 2D scattering function
+ * @param q_x: value of Q along x
+ * @param q_y: value of Q along y
+ * @return: function value
+ */
+double EllipsoidModel :: operator()(double qx, double qy) {
+	EllipsoidParameters dp;
+	// Fill parameter array
+	dp.scale      = scale();
+	dp.radius_a   = radius_a();
+	dp.radius_b   = radius_b();
+	dp.contrast   = contrast();
+	dp.background = 0.0;
+	dp.axis_theta = axis_theta();
+	dp.axis_phi   = axis_phi();
+
+	// Get the dispersion points for the radius_a
+	vector<WeightPoint> weights_rad_a;
+	radius_a.get_weights(weights_rad_a);
+
+	// Get the dispersion points for the radius_b
+	vector<WeightPoint> weights_rad_b;
+	radius_b.get_weights(weights_rad_b);
+
+	// Get angular averaging for theta
+	vector<WeightPoint> weights_theta;
+	axis_theta.get_weights(weights_theta);
+
+	// Get angular averaging for phi
+	vector<WeightPoint> weights_phi;
+	axis_phi.get_weights(weights_phi);
+
+	// Perform the computation, with all weight points
+	double sum = 0.0;
+	double norm = 0.0;
+
+	// Loop over radius weight points
+	for(int i=0; i<weights_rad_a.size(); i++) {
+		dp.radius_a = weights_rad_a[i].value;
+
+
+		// Loop over length weight points
+		for(int j=0; j<weights_rad_b.size(); j++) {
+			dp.radius_b = weights_rad_b[j].value;
+
+			// Average over theta distribution
+			for(int k=0; k<weights_theta.size(); k++) {
+				dp.axis_theta = weights_theta[k].value;
+
+				// Average over phi distribution
+				for(int l=0; l<weights_phi.size(); l++) {
+					dp.axis_phi = weights_phi[l].value;
+
+					double _ptvalue = weights_rad_a[i].weight
+						* weights_rad_b[j].weight
+						* weights_theta[k].weight
+						* weights_phi[l].weight
+						* ellipsoid_analytical_2DXY(&dp, qx, qy);
+					if (weights_theta.size()>1) {
+						_ptvalue *= sin(weights_theta[k].value);
+					}
+					sum += _ptvalue;
+
+					norm += weights_rad_a[i].weight
+						* weights_rad_b[j].weight
+						* weights_theta[k].weight
+						* weights_phi[l].weight;
+
+				}
+			}
+		}
+	}
+	// Averaging in theta needs an extra normalization
+	// factor to account for the sin(theta) term in the
+	// integration (see documentation).
+	if (weights_theta.size()>1) norm = norm / asin(1.0);
+	return sum/norm + background();
+}
+
+/**
+ * Function to evaluate 2D scattering function
+ * @param pars: parameters of the cylinder
+ * @param q: q-value
+ * @param phi: angle phi
+ * @return: function value
+ */
+double EllipsoidModel :: evaluate_rphi(double q, double phi) {
+	double qx = q*cos(phi);
+	double qy = q*sin(phi);
+	return (*this).operator()(qx, qy);
+}
Index: sansmodels/src/sans/models/c_models/ellipticalcylinder.cpp
===================================================================
--- sansmodels/src/sans/models/c_models/ellipticalcylinder.cpp	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
+++ sansmodels/src/sans/models/c_models/ellipticalcylinder.cpp	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
@@ -0,0 +1,221 @@
+/**
+	This software was developed by the University of Tennessee as part of the
+	Distributed Data Analysis of Neutron Scattering Experiments (DANSE)
+	project funded by the US National Science Foundation.
+
+	If you use DANSE applications to do scientific research that leads to
+	publication, we ask that you acknowledge the use of the software with the
+	following sentence:
+
+	"This work benefited from DANSE software developed under NSF award DMR-0520547."
+
+	copyright 2008, University of Tennessee
+ */
+
+/**
+ * Scattering model classes
+ * The classes use the IGOR library found in
+ *   sansmodels/src/libigor
+ *
+ *	TODO: refactor so that we pull in the old sansmodels.c_extensions
+ */
+
+#include <math.h>
+#include "models.hh"
+#include "parameters.hh"
+#include <stdio.h>
+using namespace std;
+
+extern "C" {
+	#include "libCylinder.h"
+	#include "elliptical_cylinder.h"
+}
+
+EllipticalCylinderModel :: EllipticalCylinderModel() {
+	scale      = Parameter(1.0);
+	r_minor    = Parameter(20.0, true);
+	r_minor.set_min(0.0);
+	r_ratio    = Parameter(1.5, true);
+	r_ratio.set_min(0.0);
+	length     = Parameter(400.0, true);
+	length.set_min(0.0);
+	contrast   = Parameter(3.e-6);
+	background = Parameter(0.0);
+	cyl_theta  = Parameter(1.57, true);
+	cyl_phi    = Parameter(0.0, true);
+	cyl_psi    = Parameter(0.0, true);
+}
+
+/**
+ * Function to evaluate 1D scattering function
+ * The NIST IGOR library is used for the actual calculation.
+ * @param q: q-value
+ * @return: function value
+ */
+double EllipticalCylinderModel :: operator()(double q) {
+	double dp[6];
+
+	dp[0] = scale();
+	dp[1] = r_minor();
+	dp[2] = r_ratio();
+	dp[3] = length();
+	dp[4] = contrast();
+	dp[5] = 0.0;
+
+	// Get the dispersion points for the r_minor
+	vector<WeightPoint> weights_rad;
+	r_minor.get_weights(weights_rad);
+
+	// Get the dispersion points for the r_ratio
+	vector<WeightPoint> weights_rat;
+	r_ratio.get_weights(weights_rat);
+
+	// Get the dispersion points for the length
+	vector<WeightPoint> weights_len;
+	length.get_weights(weights_len);
+
+	// Perform the computation, with all weight points
+	double sum = 0.0;
+	double norm = 0.0;
+
+	// Loop over r_minor weight points
+	for(int i=0; i<weights_rad.size(); i++) {
+		dp[1] = weights_rad[i].value;
+
+		// Loop over r_ratio weight points
+		for(int j=0; j<weights_rat.size(); j++) {
+			dp[2] = weights_rat[j].value;
+
+			// Loop over length weight points
+			for(int k=0; k<weights_len.size(); k++) {
+				dp[3] = weights_len[k].value;
+
+				sum += weights_rad[i].weight
+					* weights_len[k].weight
+					* weights_rat[j].weight
+					* EllipCyl20(dp, q);
+				norm += weights_rad[i].weight
+				* weights_len[k].weight
+				* weights_rat[j].weight;
+			}
+		}
+	}
+	return sum/norm + background();
+}
+
+/**
+ * Function to evaluate 2D scattering function
+ * @param q_x: value of Q along x
+ * @param q_y: value of Q along y
+ * @return: function value
+ */
+double EllipticalCylinderModel :: operator()(double qx, double qy) {
+	EllipticalCylinderParameters dp;
+	// Fill parameter array
+	dp.scale      = scale();
+	dp.r_minor    = r_minor();
+	dp.r_ratio    = r_ratio();
+	dp.length     = length();
+	dp.contrast   = contrast();
+	dp.background = 0.0;
+	dp.cyl_theta  = cyl_theta();
+	dp.cyl_phi    = cyl_phi();
+	dp.cyl_psi    = cyl_psi();
+
+	// Get the dispersion points for the r_minor
+	vector<WeightPoint> weights_rad;
+	r_minor.get_weights(weights_rad);
+
+	// Get the dispersion points for the r_ratio
+	vector<WeightPoint> weights_rat;
+	r_ratio.get_weights(weights_rat);
+
+	// Get the dispersion points for the length
+	vector<WeightPoint> weights_len;
+	length.get_weights(weights_len);
+
+	// Get angular averaging for theta
+	vector<WeightPoint> weights_theta;
+	cyl_theta.get_weights(weights_theta);
+
+	// Get angular averaging for phi
+	vector<WeightPoint> weights_phi;
+	cyl_phi.get_weights(weights_phi);
+
+	// Get angular averaging for psi
+	vector<WeightPoint> weights_psi;
+	cyl_psi.get_weights(weights_psi);
+
+	// Perform the computation, with all weight points
+	double sum = 0.0;
+	double norm = 0.0;
+
+	// Loop over minor radius weight points
+	for(int i=0; i<weights_rad.size(); i++) {
+		dp.r_minor = weights_rad[i].value;
+
+
+		// Loop over length weight points
+		for(int j=0; j<weights_len.size(); j++) {
+			dp.length = weights_len[j].value;
+
+			// Loop over r_ration weight points
+			for(int m=0; m<weights_rat.size(); m++) {
+				dp.r_ratio = weights_rat[m].value;
+
+			// Average over theta distribution
+			for(int k=0; k<weights_theta.size(); k++) {
+				dp.cyl_theta = weights_theta[k].value;
+
+				// Average over phi distribution
+				for(int l=0; l<weights_phi.size(); l++) {
+					dp.cyl_phi = weights_phi[l].value;
+
+				// Average over phi distribution
+				for(int o=0; o<weights_psi.size(); o++) {
+					dp.cyl_psi = weights_psi[o].value;
+
+					double _ptvalue = weights_rad[i].weight
+						* weights_len[j].weight
+						* weights_rat[m].weight
+						* weights_theta[k].weight
+						* weights_phi[l].weight
+						* weights_psi[o].weight
+						* elliptical_cylinder_analytical_2DXY(&dp, qx, qy);
+					if (weights_theta.size()>1) {
+						_ptvalue *= sin(weights_theta[k].value);
+					}
+					sum += _ptvalue;
+
+					norm += weights_rad[i].weight
+						* weights_len[j].weight
+						* weights_rat[m].weight
+						* weights_theta[k].weight
+						* weights_phi[l].weight
+						* weights_psi[o].weight;
+
+				}
+				}
+			}
+			}
+		}
+	}
+	// Averaging in theta needs an extra normalization
+	// factor to account for the sin(theta) term in the
+	// integration (see documentation).
+	if (weights_theta.size()>1) norm = norm / asin(1.0);
+	return sum/norm + background();
+}
+
+/**
+ * Function to evaluate 2D scattering function
+ * @param pars: parameters of the cylinder
+ * @param q: q-value
+ * @param phi: angle phi
+ * @return: function value
+ */
+double EllipticalCylinderModel :: evaluate_rphi(double q, double phi) {
+	double qx = q*cos(phi);
+	double qy = q*sin(phi);
+	return (*this).operator()(qx, qy);
+}
Index: sansmodels/src/sans/models/c_models/models.hh
===================================================================
--- sansmodels/src/sans/models/c_models/models.hh	(revision af03ddd0932915d6659404ab7824e4223ea2944f)
+++ sansmodels/src/sans/models/c_models/models.hh	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
@@ -24,5 +24,4 @@
 
 class CylinderModel{
-
 public:
 	// Model parameters
@@ -34,5 +33,4 @@
 	Parameter cyl_theta;
 	Parameter cyl_phi;
-	// TODO: replace this by an array of parameters
 
 	// Constructor
@@ -45,3 +43,105 @@
 };
 
+class SphereModel{
+public:
+	// Model parameters
+	Parameter radius;
+	Parameter scale;
+	Parameter contrast;
+	Parameter background;
+
+	// Constructor
+	SphereModel();
+
+	// Operators to get I(Q)
+	double operator()(double q);
+	double operator()(double qx, double qy);
+	double evaluate_rphi(double q, double phi);
+};
+
+class CoreShellModel{
+public:
+	// Model parameters
+	Parameter radius;
+	Parameter scale;
+	Parameter thickness;
+	Parameter core_sld;
+	Parameter shell_sld;
+	Parameter solvent_sld;
+	Parameter background;
+
+	// Constructor
+	CoreShellModel();
+
+	// Operators to get I(Q)
+	double operator()(double q);
+	double operator()(double qx, double qy);
+	double evaluate_rphi(double q, double phi);
+};
+
+class CoreShellCylinderModel{
+public:
+	// Model parameters
+	Parameter radius;
+	Parameter scale;
+	Parameter thickness;
+	Parameter length;
+	Parameter core_sld;
+	Parameter shell_sld;
+	Parameter solvent_sld;
+	Parameter background;
+	Parameter axis_theta;
+	Parameter axis_phi;
+
+	// Constructor
+	CoreShellCylinderModel();
+
+	// Operators to get I(Q)
+	double operator()(double q);
+	double operator()(double qx, double qy);
+	double evaluate_rphi(double q, double phi);
+};
+
+class EllipsoidModel{
+public:
+	// Model parameters
+	Parameter radius_a;
+	Parameter scale;
+	Parameter radius_b;
+	Parameter contrast;
+	Parameter background;
+	Parameter axis_theta;
+	Parameter axis_phi;
+
+	// Constructor
+	EllipsoidModel();
+
+	// Operators to get I(Q)
+	double operator()(double q);
+	double operator()(double qx, double qy);
+	double evaluate_rphi(double q, double phi);
+};
+
+class EllipticalCylinderModel{
+public:
+	// Model parameters
+	Parameter r_minor;
+	Parameter scale;
+	Parameter r_ratio;
+	Parameter length;
+	Parameter contrast;
+	Parameter background;
+	Parameter cyl_theta;
+	Parameter cyl_phi;
+	Parameter cyl_psi;
+
+	// Constructor
+	EllipticalCylinderModel();
+
+	// Operators to get I(Q)
+	double operator()(double q);
+	double operator()(double qx, double qy);
+	double evaluate_rphi(double q, double phi);
+};
+
 #endif
Index: sansmodels/src/sans/models/c_models/notes.txt
===================================================================
--- sansmodels/src/sans/models/c_models/notes.txt	(revision fca6936ab3704433fa2327bcd55367e8ea28de72)
+++ sansmodels/src/sans/models/c_models/notes.txt	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
@@ -11,3 +11,38 @@
 	parameters and the C++ Parameter class.
 	 
- 
+ Length = 400
+I(Q=0.001) = 450.355
+I(Q=0.001) = 472.332
+I(Q=0.001) = 4723.32
+I(Qx=0.001, Qy=0.001) = 4743.56
+I(Q=0.00447,  Phi=0.7854) = 4733
+pi=3.14159
+
+I(Qx=0.01, Qy=0.01) = 2299
+I(Q=0.0141421)         = 2326.75
+ratio 0.988074 1.01207
+
+Cyl (python):  450.35506578
+Disp(pyhton):  472.332213339
+
+--------------
+Hello cylinder
+450.35506578
+449.358988524
+446.344929175
+
+{'npts': 1, 'type': 'gaussian', 'nsigmas': 2, 'width': 0.0}
+{'npts': 1, 'type': 'flat', 'width': 0.0}
+
+ q= 0.001
+New Cyl :  450.35506578
+New Disp:  472.332213339
+
+cyl(0.01, 0.01) =96.8691
+cyl(0.0141421) =232.675
+{'type': 'array'}
+cyl(0.01, 0.01) =229.9
+Hello cylinder
+
+Old Cyl :  450.35506578
+Old Disp:  472.332213339
Index: sansmodels/src/sans/models/c_models/sphere.cpp
===================================================================
--- sansmodels/src/sans/models/c_models/sphere.cpp	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
+++ sansmodels/src/sans/models/c_models/sphere.cpp	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
@@ -0,0 +1,95 @@
+/**
+	This software was developed by the University of Tennessee as part of the
+	Distributed Data Analysis of Neutron Scattering Experiments (DANSE)
+	project funded by the US National Science Foundation.
+
+	If you use DANSE applications to do scientific research that leads to
+	publication, we ask that you acknowledge the use of the software with the
+	following sentence:
+
+	"This work benefited from DANSE software developed under NSF award DMR-0520547."
+
+	copyright 2008, University of Tennessee
+ */
+
+/**
+ * Scattering model classes
+ * The classes use the IGOR library found in
+ *   sansmodels/src/libigor
+ *
+ */
+
+#include <math.h>
+#include "models.hh"
+#include "parameters.hh"
+#include <stdio.h>
+using namespace std;
+
+extern "C" {
+	#include "libSphere.h"
+}
+
+SphereModel :: SphereModel() {
+	scale      = Parameter(1.0);
+	radius     = Parameter(20.0, true);
+	radius.set_min(0.0);
+	contrast   = Parameter(3.e-6);
+	background = Parameter(0.0);
+}
+
+/**
+ * Function to evaluate 1D scattering function
+ * The NIST IGOR library is used for the actual calculation.
+ * @param q: q-value
+ * @return: function value
+ */
+double SphereModel :: operator()(double q) {
+	double dp[4];
+
+	// Fill parameter array for IGOR library
+	// Add the background after averaging
+	dp[0] = scale();
+	dp[1] = radius();
+	dp[2] = contrast();
+	dp[3] = 0.0;
+
+	// Get the dispersion points for the radius
+	vector<WeightPoint> weights_rad;
+	radius.get_weights(weights_rad);
+
+	// Perform the computation, with all weight points
+	double sum = 0.0;
+	double norm = 0.0;
+
+	// Loop over radius weight points
+	for(int i=0; i<weights_rad.size(); i++) {
+		dp[1] = weights_rad[i].value;
+
+		sum += weights_rad[i].weight
+			* SphereForm(dp, q);
+		norm += weights_rad[i].weight;
+	}
+	return sum/norm + background();
+}
+
+/**
+ * Function to evaluate 2D scattering function
+ * @param q_x: value of Q along x
+ * @param q_y: value of Q along y
+ * @return: function value
+ */
+double SphereModel :: operator()(double qx, double qy) {
+	double q = sqrt(qx*qx + qy*qy);
+	return (*this).operator()(q);
+}
+
+/**
+ * Function to evaluate 2D scattering function
+ * @param pars: parameters of the sphere
+ * @param q: q-value
+ * @param phi: angle phi
+ * @return: function value
+ */
+double SphereModel :: evaluate_rphi(double q, double phi) {
+	return (*this).operator()(q);
+}
Index: sansmodels/src/sans/models/dispersion_models.py
===================================================================
--- sansmodels/src/sans/models/dispersion_models.py	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
+++ sansmodels/src/sans/models/dispersion_models.py	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
@@ -0,0 +1,100 @@
+"""
+    This software was developed by the University of Tennessee as part of the
+    Distributed Data Analysis of Neutron Scattering Experiments (DANSE)
+    project funded by the US National Science Foundation.
+
+    If you use DANSE applications to do scientific research that leads to
+    publication, we ask that you acknowledge the use of the software with the
+    following sentence:
+
+    "This work benefited from DANSE software developed under NSF award DMR-0520547."
+
+    copyright 2008, University of Tennessee
+"""
+
+"""
+    Class definitions for python dispersion model for 
+    model parameters. These classes are bridges to the C++
+    dispersion object.
+    
+    The ArrayDispersion class takes in numpy arrays only.
+    
+    Usage:
+    These classes can be used to set the dispersion model of a SANS model
+    parameter:
+    
+        cyl = CylinderModel()
+        cyl.set_dispersion('radius', GaussianDispersion())
+    
+    
+    After the dispersion model is set, you can access it's
+    parameter through the dispersion dictionary:
+    
+        cyl.dispersion['radius']['width'] = 5.0
+    
+    TODO: For backward compatibility, the model parameters are still kept in
+    a dictionary. The next iteration of refactoring work should involve moving
+    away from value-based parameters to object-based parameter. We want to
+    store parameters as objects so that we can unify the 'params' and 'dispersion'
+    dictionaries into a single dictionary of parameter objects that hold the
+    complete information about the parameter (units, limits, dispersion model, etc...).  
+    
+    
+"""
+import sans_extension.c_models as c_models 
+
+class DispersionModel:
+    """
+        Python bridge class for a basic dispersion model 
+        class with a constant parameter value distribution
+    """
+    def __init__(self):
+        self.cdisp = c_models.new_dispersion_model()
+        
+    def set_weights(self, values, weights):
+        """
+            Set the weights of an array dispersion
+        """
+        message = "set_weights is not available for DispersionModel.\n"
+        message += "  Solution: Use an ArrayDispersion object"
+        raise "RuntimeError", message
+        
+class GaussianDispersion(DispersionModel):
+    """
+        Python bridge class for a dispersion model based 
+        on a Gaussian distribution.
+    """
+    def __init__(self):
+        self.cdisp = c_models.new_gaussian_model()
+        
+    def set_weights(self, values, weights):
+        """
+            Set the weights of an array dispersion
+        """
+        message = "set_weights is not available for GaussiantDispersion.\n"
+        message += "  Solution: Use an ArrayDispersion object"
+        raise "RuntimeError", message
+        
+class ArrayDispersion(DispersionModel):
+    """
+        Python bridge class for a dispersion model based on arrays.
+        The user has to set a weight distribution that
+        will be used in the averaging the model parameter
+        it is applied to. 
+    """
+    def __init__(self):
+        self.cdisp = c_models.new_array_model()
+        
+    def set_weights(self, values, weights):
+        """
+            Set the weights of an array dispersion
+            Only accept numpy arrays.
+            @param values: numpy array of values
+            @param weights: numpy array of weights for each value entry
+        """
+        if len(values) != len(weights):
+            raise ValueError, "ArrayDispersion.set_weights: given arrays are of different lengths"
+        
+        c_models.set_dispersion_weights(self.cdisp, values, weights)
+        
+        
Index: sansmodels/src/sans/models/test/utest_dispersity.py
===================================================================
--- sansmodels/src/sans/models/test/utest_dispersity.py	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
+++ sansmodels/src/sans/models/test/utest_dispersity.py	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
@@ -0,0 +1,435 @@
+"""
+    Unit tests for dispersion functionality of 
+    C++ model classes
+"""
+
+import unittest, math, numpy
+
+class TestCylinder(unittest.TestCase):
+    """
+        Testing C++ Cylinder model
+    """
+    def setUp(self):
+        from sans.models.CylinderModel import CylinderModel
+        self.model= CylinderModel()
+        
+        self.model.setParam('scale', 1.0)
+        self.model.setParam('radius', 20.0)
+        self.model.setParam('length', 400.0)
+        self.model.setParam('contrast', 3.e-6)
+        self.model.setParam('background', 0.0)
+        self.model.setParam('cyl_theta', 0.0)
+        self.model.setParam('cyl_phi', 0.0)
+        
+    def test_simple(self):
+        self.assertAlmostEqual(self.model.run(0.001), 450.355, 3)
+        self.assertAlmostEqual(self.model.runXY([0.001,0.001]), 452.299, 3)
+        
+    def test_constant(self):
+        from sans.models.dispersion_models import DispersionModel
+        disp = DispersionModel()
+        self.model.setParam('scale', 10.0)
+        self.model.set_dispersion('radius', disp)
+        self.model.dispersion['radius']['width'] = 5.0
+        self.model.dispersion['radius']['npts'] = 100
+        
+        self.assertAlmostEqual(self.model.run(0.001), 4527.47250339, 3)
+        self.assertAlmostEqual(self.model.runXY([0.001, 0.001]), 4546.997777604715, 3)
+        
+    def test_gaussian(self):
+        from sans.models.dispersion_models import GaussianDispersion
+        disp = GaussianDispersion()
+        self.model.set_dispersion('radius', disp)
+        self.model.dispersion['radius']['width'] = 5.0
+        self.model.dispersion['radius']['npts'] = 100
+        self.model.setParam('scale', 10.0)
+        
+        self.assertAlmostEqual(self.model.run(0.001), 4723.32213339, 3)
+        self.assertAlmostEqual(self.model.runXY([0.001,0.001]), 4743.56, 2)
+        
+    def test_array(self):
+        """
+            Perform complete rotational average and
+            compare to 1D
+        """
+        from sans.models.dispersion_models import ArrayDispersion
+        disp_ph = ArrayDispersion()
+        disp_th = ArrayDispersion()
+        
+        values_ph = numpy.zeros(100)
+        values_th = numpy.zeros(100)
+        weights   = numpy.zeros(100)
+        for i in range(100):
+            values_ph[i]=(2.0*math.pi/99.0*i)
+            values_th[i]=(math.pi/99.0*i)
+            weights[i]=(1.0)
+        
+        disp_ph.set_weights(values_ph, weights)
+        disp_th.set_weights(values_th, weights)
+        
+        self.model.set_dispersion('cyl_theta', disp_th)
+        self.model.set_dispersion('cyl_phi', disp_ph)
+        
+        val_1d = self.model.run(math.sqrt(0.0002))
+        val_2d = self.model.runXY([0.01,0.01]) 
+        
+        self.assertTrue(math.fabs(val_1d-val_2d)/val_1d < 0.02)
+        
+class TestCoreShellCylinder(unittest.TestCase):
+    """
+        Testing C++ Cylinder model
+    """
+    def setUp(self):
+        from sans.models.CoreShellCylinderModel import CoreShellCylinderModel
+        self.model= CoreShellCylinderModel()
+        
+        self.model.setParam('scale', 1.0)
+        self.model.setParam('radius', 20.0)
+        self.model.setParam('thickness', 10.0)
+        self.model.setParam('length', 400.0)
+        self.model.setParam('core_sld', 1.e-6)
+        self.model.setParam('shell_sld', 4.e-6)
+        self.model.setParam('solvent_sld', 1.e-6)
+        self.model.setParam('background', 0.0)
+        self.model.setParam('axis_theta', 0.0)
+        self.model.setParam('axis_phi', 0.0)
+        
+    def test_simple(self):
+        """
+            Test simple 1D and 2D values
+            Numbers taken from model that passed validation, before
+            the update to C++ underlying class.
+        """
+        self.assertAlmostEqual(self.model.run(0.001), 353.55013216754583, 3)
+        self.assertAlmostEqual(self.model.runXY([0.001,0.001]), 355.25355270620543, 3)
+        
+    def test_dispersion(self):
+        """
+            Test with dispersion
+        """
+        from sans.models.DisperseModel import DisperseModel
+        disp = DisperseModel(self.model, ['radius', 'thickness', 'length'], [5, 2, 50])
+        disp.setParam('n_pts', 10)
+        self.assertAlmostEqual(disp.run(0.001), 358.44062724936009, 3)
+        self.assertAlmostEqual(disp.runXY([0.001,0.001]), 360.22673635224584, 3)
+
+    def test_new_disp(self):
+        from sans.models.dispersion_models import GaussianDispersion
+        disp_rm = GaussianDispersion()
+        self.model.set_dispersion('radius', disp_rm)
+        self.model.dispersion['radius']['width'] = 5.0
+        self.model.dispersion['radius']['npts'] = 10
+
+        disp_rr = GaussianDispersion()
+        self.model.set_dispersion('thickness', disp_rr)
+        self.model.dispersion['thickness']['width'] = 2.
+        self.model.dispersion['thickness']['npts'] = 10
+
+        disp_len = GaussianDispersion()
+        self.model.set_dispersion('length', disp_len)
+        self.model.dispersion['length']['width'] = 50.0
+        self.model.dispersion['length']['npts'] = 10
+
+        self.assertAlmostEqual(self.model.run(0.001), 358.44062724936009, 3)
+        self.assertAlmostEqual(self.model.runXY([0.001,0.001]), 360.22673635224584, 3)
+        
+
+    def test_array(self):
+        """
+            Perform complete rotational average and
+            compare to 1D
+        """
+        from sans.models.dispersion_models import ArrayDispersion
+        disp_ph = ArrayDispersion()
+        disp_th = ArrayDispersion()
+        
+        values_ph = numpy.zeros(100)
+        values_th = numpy.zeros(100)
+        weights   = numpy.zeros(100)
+        for i in range(100):
+            values_ph[i]=(2.0*math.pi/99.0*i)
+            values_th[i]=(math.pi/99.0*i)
+            weights[i]=(1.0)
+        
+        disp_ph.set_weights(values_ph, weights)
+        disp_th.set_weights(values_th, weights)
+        
+        self.model.set_dispersion('axis_theta', disp_th)
+        self.model.set_dispersion('axis_phi', disp_ph)
+        
+        val_1d = self.model.run(math.sqrt(0.0002))
+        val_2d = self.model.runXY([0.01,0.01]) 
+        
+        self.assertTrue(math.fabs(val_1d-val_2d)/val_1d < 0.02)
+        
+        
+        
+class TestCoreShell(unittest.TestCase):
+    """
+        Testing C++ Cylinder model
+    """
+    def setUp(self):
+        from sans.models.CoreShellModel import CoreShellModel
+        self.model= CoreShellModel()
+        
+        self.model.setParam('scale', 1.0)
+        self.model.setParam('radius', 60.0)
+        self.model.setParam('thickness', 10.0)
+        self.model.setParam('core_sld', 1.e-6)
+        self.model.setParam('shell_sld', 2.e-6)
+        self.model.setParam('solvent_sld', 3.e-6)
+        self.model.setParam('background', 0.0)
+        
+    def test_simple(self):
+        """
+            Test simple 1D and 2D values
+            Numbers taken from model that passed validation, before
+            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)
+        
+    def test_dispersion(self):
+        """
+            Test with dispersion
+        """
+        from sans.models.DisperseModel import DisperseModel
+        disp = DisperseModel(self.model, ['radius', 'thickness'], [10, 2])
+        disp.setParam('n_pts', 10)
+        self.assertAlmostEqual(disp.run(0.001), 407.344127907553, 3)
+    
+    def test_new_disp(self):
+        from sans.models.dispersion_models import GaussianDispersion
+        disp_rm = GaussianDispersion()
+        self.model.set_dispersion('radius', disp_rm)
+        self.model.dispersion['radius']['width'] = 10.
+        self.model.dispersion['radius']['npts'] = 10
+
+        disp_rr = GaussianDispersion()
+        self.model.set_dispersion('thickness', disp_rr)
+        self.model.dispersion['thickness']['width'] = 2.
+        self.model.dispersion['thickness']['npts'] = 10
+
+
+        self.assertAlmostEqual(self.model.run(0.001), 407.344127907553, 3)
+
+        
+class TestEllipsoid(unittest.TestCase):
+    """
+        Testing C++ Cylinder model
+    """
+    def setUp(self):
+        from sans.models.EllipsoidModel import EllipsoidModel
+        self.model= EllipsoidModel()
+        
+        self.model.setParam('scale', 1.0)
+        self.model.setParam('radius_a', 20.0)
+        self.model.setParam('radius_b', 400.0)
+        self.model.setParam('contrast', 3.e-6)
+        self.model.setParam('background', 0.0)
+        self.model.setParam('axis_theta', 1.57)
+        self.model.setParam('axis_phi', 0.0)
+        
+    def test_simple(self):
+        """
+            Test simple 1D and 2D values
+            Numbers taken from model that passed validation, before
+            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)
+
+    def test_dispersion(self):
+        """
+            Test with dispersion
+        """
+        from sans.models.DisperseModel import DisperseModel
+        disp = DisperseModel(self.model, ['radius_a', 'radius_b'], [5, 50])
+        disp.setParam('n_pts', 10)
+        self.assertAlmostEqual(disp.run(0.001), 11948.72581312305, 3)
+        self.assertAlmostEqual(disp.runXY([0.001,0.001]), 11811.972359807551, 3)
+        
+    def test_new_disp(self):
+        from sans.models.dispersion_models import GaussianDispersion
+        disp_rm = GaussianDispersion()
+        self.model.set_dispersion('radius_a', disp_rm)
+        self.model.dispersion['radius_a']['width'] = 5.0
+        self.model.dispersion['radius_a']['npts'] = 10
+
+        disp_rr = GaussianDispersion()
+        self.model.set_dispersion('radius_b', disp_rr)
+        self.model.dispersion['radius_b']['width'] = 50.
+        self.model.dispersion['radius_b']['npts'] = 10
+
+
+        self.assertAlmostEqual(self.model.run(0.001), 11948.72581312305, 3)
+        self.assertAlmostEqual(self.model.runXY([0.001,0.001]), 11811.972359807551, 3)
+
+    def test_array(self):
+        """
+            Perform complete rotational average and
+            compare to 1D
+        """
+        from sans.models.dispersion_models import ArrayDispersion
+        disp_ph = ArrayDispersion()
+        disp_th = ArrayDispersion()
+        
+        values_ph = numpy.zeros(100)
+        values_th = numpy.zeros(100)
+        weights   = numpy.zeros(100)
+        for i in range(100):
+            values_ph[i]=(2.0*math.pi/99.0*i)
+            values_th[i]=(math.pi/99.0*i)
+            weights[i]=(1.0)
+        
+        disp_ph.set_weights(values_ph, weights)
+        disp_th.set_weights(values_th, weights)
+        
+        self.model.set_dispersion('axis_theta', disp_th)
+        self.model.set_dispersion('axis_phi', disp_ph)
+        
+        val_1d = self.model.run(math.sqrt(0.0002))
+        val_2d = self.model.runXY([0.01,0.01]) 
+        
+        self.assertTrue(math.fabs(val_1d-val_2d)/val_1d < 0.02)
+        
+        
+        
+class TestSphere(unittest.TestCase):
+    """
+        Testing C++ Cylinder model
+    """
+    def setUp(self):
+        from sans.models.SphereModel import SphereModel
+        self.model= SphereModel()
+        
+        self.model.setParam('scale', 1.0)
+        self.model.setParam('radius', 60.0)
+        self.model.setParam('contrast', 1)
+        self.model.setParam('background', 0.0)
+        
+    def test_simple(self):
+        """
+            Test simple 1D and 2D values
+            Numbers taken from model that passed validation, before
+            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, 3)
+
+    def test_dispersion(self):
+        """
+            Test with dispersion
+        """
+        from sans.models.DisperseModel import DisperseModel
+        disp = DisperseModel(self.model, ['radius'], [10])
+        disp.setParam('n_pts', 10)
+        self.assertAlmostEqual(disp.run(0.001), 96795008379475.25, 3)
+        
+    def test_new_disp(self):
+        from sans.models.dispersion_models import GaussianDispersion
+        disp_rm = GaussianDispersion()
+        self.model.set_dispersion('radius', disp_rm)
+        self.model.dispersion['radius']['width'] = 10.0
+        self.model.dispersion['radius']['npts'] = 10
+
+        self.assertAlmostEqual(self.model.run(0.001), 96795008379475.25, 3)
+        
+        
+        
+class TestEllipticalCylinder(unittest.TestCase):
+    """
+        Testing C++ Cylinder model
+    """
+    def setUp(self):
+        from sans.models.EllipticalCylinderModel import EllipticalCylinderModel
+        self.model= EllipticalCylinderModel()
+        
+        self.model.setParam('scale', 1.0)
+        self.model.setParam('r_minor', 20.0)
+        self.model.setParam('r_ratio', 1.5)
+        self.model.setParam('length', 400.0)
+        self.model.setParam('contrast', 3.e-6)
+        self.model.setParam('background', 0.0)
+        self.model.setParam('cyl_theta', 1.57)
+        self.model.setParam('cyl_phi', 0.0)
+        self.model.setParam('cyl_psi', 0.0)
+        
+    def test_simple(self):
+        """
+            Test simple 1D and 2D values
+            Numbers taken from model that passed validation, before
+            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, 3)
+        
+    def test_dispersion(self):
+        """
+            Test with dispersion
+        """
+        from sans.models.DisperseModel import DisperseModel
+        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)
+        self.assertAlmostEqual(disp.runXY([0.001,0.001]), 704.63525988095705, 3)
+
+    def test_new_disp(self):
+        from sans.models.dispersion_models import GaussianDispersion
+        disp_rm = GaussianDispersion()
+        self.model.set_dispersion('r_minor', disp_rm)
+        self.model.dispersion['r_minor']['width'] = 5.0
+        self.model.dispersion['r_minor']['npts'] = 10
+
+        disp_rr = GaussianDispersion()
+        self.model.set_dispersion('r_ratio', disp_rr)
+        self.model.dispersion['r_ratio']['width'] = 0.25
+        self.model.dispersion['r_ratio']['npts'] = 10
+
+        disp_len = GaussianDispersion()
+        self.model.set_dispersion('length', disp_len)
+        self.model.dispersion['length']['width'] = 50.0
+        self.model.dispersion['length']['npts'] = 10
+
+        self.assertAlmostEqual(self.model.run(0.001), 711.18048194151925, 3)
+        self.assertAlmostEqual(self.model.runXY([0.001,0.001]), 704.63525988095705, 3)
+        
+
+    def test_array(self):
+        """
+            Perform complete rotational average and
+            compare to 1D
+        """
+        from sans.models.dispersion_models import ArrayDispersion
+        disp_ph = ArrayDispersion()
+        disp_th = ArrayDispersion()
+        disp_ps = ArrayDispersion()
+        
+        values_ph = numpy.zeros(100)
+        values_th = numpy.zeros(100)
+        values_ps = numpy.zeros(100)
+        weights   = numpy.zeros(100)
+        for i in range(100):
+            values_ps[i]=(2.0*math.pi/99.0*i)
+            values_ph[i]=(2.0*math.pi/99.0*i)
+            values_th[i]=(math.pi/99.0*i)
+            weights[i]=(1.0)
+        
+        disp_ph.set_weights(values_ph, weights)
+        disp_th.set_weights(values_th, weights)
+        disp_ps.set_weights(values_ps, weights)
+        
+        self.model.set_dispersion('cyl_theta', disp_th)
+        self.model.set_dispersion('cyl_phi', disp_ph)
+        self.model.set_dispersion('cyl_psi', disp_ps)
+        
+        val_1d = self.model.run(math.sqrt(0.0002))
+        val_2d = self.model.runXY([0.01,0.01]) 
+        
+        self.assertTrue(math.fabs(val_1d-val_2d)/val_1d < 0.02)
+        
+        
+  
+if __name__ == '__main__':
+    unittest.main()
+   
Index: sansmodels/src/setup.py
===================================================================
--- sansmodels/src/setup.py	(revision af03ddd0932915d6659404ab7824e4223ea2944f)
+++ sansmodels/src/setup.py	(revision 0f5bc9ff2027432913fd8d2cb067675244acdaf6)
@@ -93,6 +93,8 @@
      sources = [
         srcdir+"/c_models.cpp",
-        srcdir+"/CSphereModel.c",
-        srcdir+"/sphere.c",
+        #srcdir+"/CSphereModel.c",
+        #srcdir+"/sphere.c",
+        "sans/models/c_models/CSphereModel.cpp",
+        "sans/models/c_models/sphere.cpp",
         #srcdir+"/CCylinderModel.c",
         "sans/models/c_models/CCylinderModel.cpp",
@@ -101,11 +103,19 @@
         "sans/models/c_models/dispersion_visitor.cpp",
         srcdir+"/cylinder.c",
-        srcdir+"/CCoreShellCylinderModel.c",
+        #srcdir+"/CCoreShellCylinderModel.c",
+        "sans/models/c_models/CCoreShellCylinderModel.cpp",
+        "sans/models/c_models/coreshellcylinder.cpp",
         srcdir+"/core_shell_cylinder.c",
-        srcdir+"/CCoreShellModel.c",
-        srcdir+"/core_shell.c",
-        srcdir+"/CEllipsoidModel.c",
+        #srcdir+"/CCoreShellModel.c",
+        #srcdir+"/core_shell.c",
+        "sans/models/c_models/CCoreShellModel.cpp",
+        "sans/models/c_models/coreshellsphere.cpp",
+        #srcdir+"/CEllipsoidModel.c",
+        "sans/models/c_models/CEllipsoidModel.cpp",
+        "sans/models/c_models/ellipsoid.cpp",        
         srcdir+"/ellipsoid.c",
-        srcdir+"/CEllipticalCylinderModel.c",
+        #srcdir+"/CEllipticalCylinderModel.c",
+        "sans/models/c_models/CEllipticalCylinderModel.cpp",
+        "sans/models/c_models/ellipticalcylinder.cpp",                
         srcdir+"/elliptical_cylinder.c",
         srcdir+"/disperser.c",
