Index: sansmodels/src/c_models/corefourshell.cpp
===================================================================
--- sansmodels/src/c_models/corefourshell.cpp	(revision 67424cdc019bde2f4cbc089828ff660e6f914c2c)
+++ sansmodels/src/c_models/corefourshell.cpp	(revision 011e0e4569d018a5f39b70032285911f67750383)
@@ -21,13 +21,28 @@
 
 #include <math.h>
-#include "models.hh"
 #include "parameters.hh"
 #include <stdio.h>
 using namespace std;
+#include "corefourshell.h"
 
 extern "C" {
 	#include "libSphere.h"
-	#include "corefourshell.h"
-}
+}
+
+typedef struct {
+  double scale;
+  double rad_core0;
+  double sld_core0;
+  double thick_shell1;
+  double sld_shell1;
+  double thick_shell2;
+  double sld_shell2;
+  double thick_shell3;
+  double sld_shell3;
+  double thick_shell4;
+  double sld_shell4;
+  double sld_solv;
+  double background;
+} CoreFourShellParameters;
 
 CoreFourShellModel :: CoreFourShellModel() {
Index: sansmodels/src/c_models/coreshellcylinder.cpp
===================================================================
--- sansmodels/src/c_models/coreshellcylinder.cpp	(revision 67424cdc019bde2f4cbc089828ff660e6f914c2c)
+++ sansmodels/src/c_models/coreshellcylinder.cpp	(revision 011e0e4569d018a5f39b70032285911f67750383)
@@ -18,18 +18,104 @@
  *   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;
+#include "core_shell_cylinder.h"
 
 extern "C" {
 	#include "libCylinder.h"
 	#include "libStructureFactor.h"
-	#include "core_shell_cylinder.h"
-}
+}
+
+typedef struct {
+  double scale;
+  double radius;
+  double thickness;
+  double length;
+  double core_sld;
+  double shell_sld;
+  double solvent_sld;
+  double background;
+  double axis_theta;
+  double axis_phi;
+} CoreShellCylinderParameters;
+
+
+/**
+ * Function to evaluate 2D scattering function
+ * @param pars: parameters of the core-shell cylinder
+ * @param q: q-value
+ * @param q_x: q_x / q
+ * @param q_y: q_y / q
+ * @return: function value
+ */
+static double core_shell_cylinder_analytical_2D_scaled(CoreShellCylinderParameters *pars, double q, double q_x, double q_y) {
+  double cyl_x, cyl_y, cyl_z;
+  double q_z;
+  double alpha, vol, cos_val;
+  double answer;
+  //convert angle degree to radian
+  double pi = 4.0*atan(1.0);
+  double theta = pars->axis_theta * pi/180.0;
+  double phi = pars->axis_phi * pi/180.0;
+
+    // Cylinder orientation
+    cyl_x = sin(theta) * cos(phi);
+    cyl_y = sin(theta) * sin(phi);
+    cyl_z = cos(theta);
+
+    // q vector
+    q_z = 0;
+
+    // Compute the angle btw vector q and the
+    // axis of the cylinder
+    cos_val = cyl_x*q_x + cyl_y*q_y + cyl_z*q_z;
+
+    // The following test should always pass
+    if (fabs(cos_val)>1.0) {
+      printf("core_shell_cylinder_analytical_2D: Unexpected error: cos(alpha)=%g\n", cos_val);
+      return 0;
+    }
+
+  alpha = acos( cos_val );
+
+  // Call the IGOR library function to get the kernel
+  answer = CoreShellCylKernel(q, pars->radius, pars->thickness,
+                pars->core_sld,pars->shell_sld,
+                pars->solvent_sld, pars->length/2.0, alpha) / fabs(sin(alpha));
+
+  //normalize by cylinder volume
+  vol=pi*(pars->radius+pars->thickness)
+      *(pars->radius+pars->thickness)
+      *(pars->length+2.0*pars->thickness);
+  answer /= vol;
+
+  //convert to [cm-1]
+  answer *= 1.0e8;
+
+  //Scale
+  answer *= pars->scale;
+
+  // add in the background
+  answer += pars->background;
+
+  return answer;
+}
+
+/**
+ * Function to evaluate 2D scattering function
+ * @param pars: parameters of the core-shell cylinder
+ * @param q: q-value
+ * @return: function value
+ */
+static double core_shell_cylinder_analytical_2DXY(CoreShellCylinderParameters *pars, double qx, double qy) {
+  double q;
+  q = sqrt(qx*qx+qy*qy);
+    return core_shell_cylinder_analytical_2D_scaled(pars, q, qx/q, qy/q);
+}
+
 
 CoreShellCylinderModel :: CoreShellCylinderModel() {
Index: sansmodels/src/c_models/coreshellsphere.cpp
===================================================================
--- sansmodels/src/c_models/coreshellsphere.cpp	(revision 67424cdc019bde2f4cbc089828ff660e6f914c2c)
+++ sansmodels/src/c_models/coreshellsphere.cpp	(revision 011e0e4569d018a5f39b70032285911f67750383)
@@ -21,24 +21,33 @@
 
 #include <math.h>
-#include "models.hh"
 #include "parameters.hh"
 #include <stdio.h>
 using namespace std;
+#include "core_shell.h"
 
 extern "C" {
-	#include "libSphere.h"
-	#include "core_shell.h"
+#include "libSphere.h"
 }
 
+typedef struct {
+  double scale;
+  double radius;
+  double thickness;
+  double core_sld;
+  double shell_sld;
+  double solvent_sld;
+  double background;
+} CoreShellParameters;
+
 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);
+  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);
 }
 
@@ -50,55 +59,55 @@
  */
 double CoreShellModel :: operator()(double q) {
-	double dp[7];
+  double dp[7];
 
-	// Fill parameter array for IGOR library
-	// Add the background after averaging
+  // 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;
+  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 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 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;
-	double vol = 0.0;
+  // Perform the computation, with all weight points
+  double sum = 0.0;
+  double norm = 0.0;
+  double vol = 0.0;
 
-	// Loop over radius weight points
-	for(size_t i=0; i<weights_rad.size(); i++) {
-		dp[1] = weights_rad[i].value;
+  // Loop over radius weight points
+  for(size_t i=0; i<weights_rad.size(); i++) {
+    dp[1] = weights_rad[i].value;
 
-		// Loop over thickness weight points
-		for(size_t j=0; j<weights_thick.size(); j++) {
-			dp[2] = weights_thick[j].value;
-			//Un-normalize SphereForm by volume
-			sum += weights_rad[i].weight
-				* weights_thick[j].weight * CoreShellForm(dp, q)* pow(weights_rad[i].value+weights_thick[j].value,3);
+    // Loop over thickness weight points
+    for(size_t j=0; j<weights_thick.size(); j++) {
+      dp[2] = weights_thick[j].value;
+      //Un-normalize SphereForm by volume
+      sum += weights_rad[i].weight
+          * weights_thick[j].weight * CoreShellForm(dp, q)* pow(weights_rad[i].value+weights_thick[j].value,3);
 
-			//Find average volume
-			vol += weights_rad[i].weight * weights_thick[j].weight
-				* pow(weights_rad[i].value+weights_thick[j].value,3);
-			norm += weights_rad[i].weight
-				* weights_thick[j].weight;
-		}
-	}
+      //Find average volume
+      vol += weights_rad[i].weight * weights_thick[j].weight
+          * pow(weights_rad[i].value+weights_thick[j].value,3);
+      norm += weights_rad[i].weight
+          * weights_thick[j].weight;
+    }
+  }
 
-	if (vol != 0.0 && norm != 0.0) {
-		//Re-normalize by avg volume
-		sum = sum/(vol/norm);}
+  if (vol != 0.0 && norm != 0.0) {
+    //Re-normalize by avg volume
+    sum = sum/(vol/norm);}
 
-	return sum/norm + background();
+  return sum/norm + background();
 }
 
@@ -110,6 +119,6 @@
  */
 double CoreShellModel :: operator()(double qx, double qy) {
-	double q = sqrt(qx*qx + qy*qy);
-	return (*this).operator()(q);
+  double q = sqrt(qx*qx + qy*qy);
+  return (*this).operator()(q);
 }
 
@@ -122,5 +131,5 @@
  */
 double CoreShellModel :: evaluate_rphi(double q, double phi) {
-	return (*this).operator()(q);
+  return (*this).operator()(q);
 }
 /**
@@ -129,41 +138,41 @@
  */
 double CoreShellModel :: calculate_ER() {
-	CoreShellParameters dp;
+  CoreShellParameters dp;
 
-	dp.radius     = radius();
-	dp.thickness  = thickness();
+  dp.radius     = radius();
+  dp.thickness  = thickness();
 
-	double rad_out = 0.0;
+  double rad_out = 0.0;
 
-	// Perform the computation, with all weight points
-	double sum = 0.0;
-	double norm = 0.0;
+  // Perform the computation, with all weight points
+  double sum = 0.0;
+  double norm = 0.0;
 
 
-	// Get the dispersion points for the major shell
-	vector<WeightPoint> weights_thickness;
-	thickness.get_weights(weights_thickness);
+  // Get the dispersion points for the major shell
+  vector<WeightPoint> weights_thickness;
+  thickness.get_weights(weights_thickness);
 
-	// Get the dispersion points for the minor shell
-	vector<WeightPoint> weights_radius ;
-	radius.get_weights(weights_radius);
+  // Get the dispersion points for the minor shell
+  vector<WeightPoint> weights_radius ;
+  radius.get_weights(weights_radius);
 
-	// Loop over major shell weight points
-	for(int j=0; j< (int)weights_thickness.size(); j++) {
-		dp.thickness = weights_thickness[j].value;
-		for(int k=0; k< (int)weights_radius.size(); k++) {
-			dp.radius = weights_radius[k].value;
-			sum += weights_thickness[j].weight
-				* weights_radius[k].weight*(dp.radius+dp.thickness);
-			norm += weights_thickness[j].weight* weights_radius[k].weight;
-		}
-	}
-	if (norm != 0){
-		//return the averaged value
-		rad_out =  sum/norm;}
-	else{
-		//return normal value
-		rad_out = (dp.radius+dp.thickness);}
+  // Loop over major shell weight points
+  for(int j=0; j< (int)weights_thickness.size(); j++) {
+    dp.thickness = weights_thickness[j].value;
+    for(int k=0; k< (int)weights_radius.size(); k++) {
+      dp.radius = weights_radius[k].value;
+      sum += weights_thickness[j].weight
+          * weights_radius[k].weight*(dp.radius+dp.thickness);
+      norm += weights_thickness[j].weight* weights_radius[k].weight;
+    }
+  }
+  if (norm != 0){
+    //return the averaged value
+    rad_out =  sum/norm;}
+  else{
+    //return normal value
+    rad_out = (dp.radius+dp.thickness);}
 
-	return rad_out;
+  return rad_out;
 }
Index: sansmodels/src/c_models/ellipsoid.cpp
===================================================================
--- sansmodels/src/c_models/ellipsoid.cpp	(revision 67424cdc019bde2f4cbc089828ff660e6f914c2c)
+++ sansmodels/src/c_models/ellipsoid.cpp	(revision 011e0e4569d018a5f39b70032285911f67750383)
@@ -18,17 +18,99 @@
  *   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;
+#include "ellipsoid.h"
 
 extern "C" {
 	#include "libCylinder.h"
 	#include "libStructureFactor.h"
-	#include "ellipsoid.h"
+}
+
+typedef struct {
+  double scale;
+  double radius_a;
+  double radius_b;
+  double sldEll;
+  double sldSolv;
+  double background;
+  double axis_theta;
+  double axis_phi;
+} EllipsoidParameters;
+
+/**
+ * Function to evaluate 2D scattering function
+ * @param pars: parameters of the ellipsoid
+ * @param q: q-value
+ * @param q_x: q_x / q
+ * @param q_y: q_y / q
+ * @return: function value
+ */
+static double ellipsoid_analytical_2D_scaled(EllipsoidParameters *pars, double q, double q_x, double q_y) {
+  double cyl_x, cyl_y, cyl_z;
+  double q_z;
+  double alpha, vol, cos_val;
+  double answer;
+  //convert angle degree to radian
+  double pi = 4.0*atan(1.0);
+  double theta = pars->axis_theta * pi/180.0;
+  double phi = pars->axis_phi * pi/180.0;
+
+    // Ellipsoid orientation
+    cyl_x = sin(theta) * cos(phi);
+    cyl_y = sin(theta) * sin(phi);
+    cyl_z = cos(theta);
+
+    // q vector
+    q_z = 0.0;
+
+    // Compute the angle btw vector q and the
+    // axis of the cylinder
+    cos_val = cyl_x*q_x + cyl_y*q_y + cyl_z*q_z;
+
+    // The following test should always pass
+    if (fabs(cos_val)>1.0) {
+      printf("ellipsoid_ana_2D: Unexpected error: cos(alpha)>1\n");
+      return 0;
+    }
+
+    // Angle between rotation axis and q vector
+  alpha = acos( cos_val );
+
+  // Call the IGOR library function to get the kernel
+  answer = EllipsoidKernel(q, pars->radius_b, pars->radius_a, cos_val);
+
+  // Multiply by contrast^2
+  answer *= (pars->sldEll - pars->sldSolv) * (pars->sldEll - pars->sldSolv);
+
+  //normalize by cylinder volume
+    vol = 4.0/3.0 * acos(-1.0) * pars->radius_b * pars->radius_b * pars->radius_a;
+  answer *= vol;
+
+  //convert to [cm-1]
+  answer *= 1.0e8;
+
+  //Scale
+  answer *= pars->scale;
+
+  // add in the background
+  answer += pars->background;
+
+  return answer;
+}
+
+/**
+ * Function to evaluate 2D scattering function
+ * @param pars: parameters of the ellipsoid
+ * @param q: q-value
+ * @return: function value
+ */
+static double ellipsoid_analytical_2DXY(EllipsoidParameters *pars, double qx, double qy) {
+  double q;
+  q = sqrt(qx*qx+qy*qy);
+    return ellipsoid_analytical_2D_scaled(pars, q, qx/q, qy/q);
 }
 
Index: sansmodels/src/c_models/models.hh
===================================================================
--- sansmodels/src/c_models/models.hh	(revision 0ba3b08c35684feff124c0657dad569d13717a18)
+++ sansmodels/src/c_models/models.hh	(revision 011e0e4569d018a5f39b70032285911f67750383)
@@ -27,100 +27,4 @@
 using namespace std;
 
-
-
-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 calculate_ER();
-	double evaluate_rphi(double q, double phi);
-};
-
-class CoreFourShellModel{
-public:
-	// Model parameters
-	Parameter scale;
-	Parameter rad_core0;
-	Parameter sld_core0;
-	Parameter thick_shell1;
-	Parameter sld_shell1;
-	Parameter thick_shell2;
-	Parameter sld_shell2;
-	Parameter thick_shell3;
-	Parameter sld_shell3;
-	Parameter thick_shell4;
-	Parameter sld_shell4;
-	Parameter sld_solv;
-	Parameter background;
-
-	// Constructor
-	CoreFourShellModel();
-
-	// Operators to get I(Q)
-	double operator()(double q);
-	double operator()(double qx, double qy);
-	double calculate_ER();
-	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 calculate_ER();
-	double evaluate_rphi(double q, double phi);
-};
-
-class EllipsoidModel{
-public:
-	// Model parameters
-	Parameter radius_a;
-	Parameter scale;
-	Parameter radius_b;
-	Parameter sldEll;
-	Parameter sldSolv;
-	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 calculate_ER();
-	double evaluate_rphi(double q, double phi);
-};
-
 class EllipticalCylinderModel{
 public:
