Index: sasmodels/models/core_shell_bicelle.c
===================================================================
--- sasmodels/models/core_shell_bicelle.c	(revision 8007311db20a23588c62551c8a69f089cf8835e1)
+++ sasmodels/models/core_shell_bicelle.c	(revision 8007311db20a23588c62551c8a69f089cf8835e1)
@@ -0,0 +1,209 @@
+double form_volume(double radius, double rim_thickness, double face_thickness, double length);
+double Iq(double q,
+          double radius,
+          double rim_thickness,
+          double face_thickness,
+          double length,
+          double core_sld,
+          double face_sld,
+          double rim_sld,
+          double solvent_sld);
+
+
+double Iqxy(double qx, double qy,
+          double radius,
+          double rim_thickness,
+          double face_thickness,
+          double length,
+          double core_sld,
+          double face_sld,
+          double rim_sld,
+          double solvent_sld,
+          double theta,
+          double phi);
+
+
+double form_volume(double radius, double rim_thickness, double face_thickness, double length)
+{
+    return M_PI*(radius+rim_thickness)*(radius+rim_thickness)*(length+2*face_thickness);
+}
+
+static double
+bicelle_kernel(double qq,
+              double rad,
+              double radthick,
+              double facthick,
+              double length,
+              double rhoc,
+              double rhoh,
+              double rhor,
+              double rhosolv,
+              double dum)
+{
+	double dr1,dr2,dr3;
+	double besarg1,besarg2;
+	double vol1,vol2,vol3;
+	double sinarg1,sinarg2;
+	double t1,t2,t3;
+	double retval,si1,si2,be1,be2;
+
+	dr1 = rhoc-rhoh;
+	dr2 = rhor-rhosolv;
+	dr3=  rhoh-rhor;
+	vol1 = M_PI*rad*rad*(2.0*length);
+	vol2 = M_PI*(rad+radthick)*(rad+radthick)*(2.0*length+2.0*facthick);
+	vol3= M_PI*(rad)*(rad)*(2.0*length+2.0*facthick);
+	besarg1 = qq*rad*sin(dum);
+	besarg2 = qq*(rad+radthick)*sin(dum);
+	sinarg1 = qq*length*cos(dum);
+	sinarg2 = qq*(length+facthick)*cos(dum);
+
+	if(besarg1 == 0) {
+		be1 = 0.5;
+	} else {
+		be1 = J1(besarg1)/besarg1;
+	}
+	if(besarg2 == 0) {
+		be2 = 0.5;
+	} else {
+		be2 = J1(besarg2)/besarg2;
+	}
+	if(sinarg1 == 0) {
+		si1 = 1.0;
+	} else {
+		si1 = sin(sinarg1)/sinarg1;
+	}
+	if(sinarg2 == 0) {
+		si2 = 1.0;
+	} else {
+		si2 = sin(sinarg2)/sinarg2;
+	}
+	t1 = 2.0*vol1*dr1*si1*be1;
+	t2 = 2.0*vol2*dr2*si2*be2;
+	t3 = 2.0*vol3*dr3*si2*be1;
+
+	retval = ((t1+t2+t3)*(t1+t2+t3))*sin(dum);
+	return(retval);
+
+}
+
+static double
+bicelle_integration(double qq,
+                   double rad,
+                   double radthick,
+                   double facthick,
+                   double length,
+                   double rhoc,
+                   double rhoh,
+                   double rhor,
+                   double rhosolv)
+{
+
+
+	double answer,halfheight;
+	double lolim,uplim,summ,yyy,zi;
+	int nord,i;
+
+	// set up the integration end points
+	nord = 76;
+	lolim = 0.0;
+	uplim = M_PI/2;
+	halfheight = length/2.0;
+
+	summ = 0.0;
+	i=0;
+	for(i=0;i<nord;i++) {
+		zi = ( Gauss76Z[i]*(uplim-lolim) + uplim + lolim )/2.0;
+		yyy = Gauss76Wt[i] * bicelle_kernel(qq, rad, radthick, facthick,
+		                     halfheight, rhoc, rhoh, rhor,rhosolv, zi);
+		summ += yyy;
+	}
+
+	// calculate value of integral to return
+	answer = (uplim-lolim)/2.0*summ;
+	return(answer);
+}
+
+static double
+bicelle_kernel_2d(double q, double q_x, double q_y,
+          double radius,
+          double rim_thickness,
+          double face_thickness,
+          double length,
+          double core_sld,
+          double face_sld,
+          double rim_sld,
+          double solvent_sld,
+          double theta,
+          double phi)
+{
+    double cyl_x, cyl_y;
+    double alpha, cos_val;
+    double answer;
+
+    //convert angle degree to radian
+    theta *= M_PI/180.0;
+    phi *= M_PI/180.0;
+
+    // Cylinder orientation
+    cyl_x = cos(theta) * cos(phi);
+    cyl_y = sin(theta);
+
+    // Compute the angle btw vector q and the axis of the cylinder
+    cos_val = cyl_x*q_x + cyl_y*q_y;
+    alpha = acos( cos_val );
+
+    // Get the kernel
+    answer = bicelle_kernel(q, radius, rim_thickness, face_thickness,
+                           length/2.0, core_sld, face_sld, rim_sld,
+                           solvent_sld, alpha) / fabs(sin(alpha));
+
+    answer *= 1.0e-4;
+
+    return answer;
+}
+
+double Iq(double q,
+          double radius,
+          double rim_thickness,
+          double face_thickness,
+          double length,
+          double core_sld,
+          double face_sld,
+          double rim_sld,
+          double solvent_sld)
+{
+    double intensity = bicelle_integration(q, radius, rim_thickness, face_thickness,
+                       length, core_sld, face_sld, rim_sld, solvent_sld);
+    return intensity*1.0e-4;
+}
+
+
+double Iqxy(double qx, double qy,
+          double radius,
+          double rim_thickness,
+          double face_thickness,
+          double length,
+          double core_sld,
+          double face_sld,
+          double rim_sld,
+          double solvent_sld,
+          double theta,
+          double phi)
+{
+    double q;
+    q = sqrt(qx*qx+qy*qy);
+    double intensity = bicelle_kernel_2d(q, qx/q, qy/q,
+                      radius,
+                      rim_thickness,
+                      face_thickness,
+                      length,
+                      core_sld,
+                      face_sld,
+                      rim_sld,
+                      solvent_sld,
+                      theta,
+                      phi);
+
+    return intensity;
+}
Index: sasmodels/models/core_shell_bicelle.py
===================================================================
--- sasmodels/models/core_shell_bicelle.py	(revision 8007311db20a23588c62551c8a69f089cf8835e1)
+++ sasmodels/models/core_shell_bicelle.py	(revision 8007311db20a23588c62551c8a69f089cf8835e1)
@@ -0,0 +1,151 @@
+r"""
+
+Definition
+----------
+This model provides the form factor for a circular cylinder with a core-shell
+scattering length density profile.
+The form factor is normalized by the particle volume.
+
+.. _core-shell-cylinder-geometry:
+
+.. figure:: img/core_shell_bicelle_geometry.png
+
+    (Graphic from DOI: 10.1039/C0NP00002G)
+
+
+The output of the 1D scattering intensity function for randomly oriented
+cylinders is then given by the equation above.
+
+The *theta* and *phi* parameters are not used for the 1D output.
+Our implementation of the scattering kernel and the 1D scattering intensity
+use the c-library from NIST.
+
+.. figure:: img/core_shell_bicelle_1d.jpg
+
+    1D plot using the default values (w/200 data point).
+
+.. figure:: img/core_shell_bicelle_fig1.jpg
+
+    Definition of the angles for the oriented CoreShellBicelleModel.
+
+.. figure:: img/core_shell_bicelle_fig2.jpg
+
+    Examples of the angles for oriented pp against the detector plane.
+
+References
+----------
+
+L A Feigin and D I Svergun,
+*Structure Analysis by Small-Angle X-Ray and Neutron Scattering,*
+Plenum Press, New York, (1987)
+
+"""
+
+from numpy import inf, sin, cos
+
+name = "core_shell_bicelle"
+title = "Circular cylinder with a core-shell scattering length density profile.."
+description = """
+    P(q,alpha)= scale/Vs*f(q)^(2) + bkg,  where: f(q)= 2(core_sld
+    - solvant_sld)* Vc*sin[qLcos(alpha/2)]
+    /[qLcos(alpha/2)]*J1(qRsin(alpha))
+    /[qRsin(alpha)]+2(shell_sld-solvent_sld)
+    *Vs*sin[q(L+T)cos(alpha/2)][[q(L+T)
+    *cos(alpha/2)]*J1(q(R+T)sin(alpha))
+    /q(R+T)sin(alpha)]
+
+    alpha:is the angle between the axis of
+    the cylinder and the q-vector
+    Vs: the volume of the outer shell
+    Vc: the volume of the core
+    L: the length of the core
+    shell_sld: the scattering length density
+    of the shell
+    solvent_sld: the scattering length density
+    of the solvent
+    bkg: the background
+    T: the thickness
+    R+T: is the outer radius
+    L+2T: The total length of the outershell
+    J1: the first order Bessel function
+    theta: axis_theta of the cylinder
+    phi: the axis_phi of the cylinder...
+        """
+category = "shape:cylinder"
+
+# pylint: disable=bad-whitespace, line-too-long
+#             ["name", "units", default, [lower, upper], "type", "description"],
+parameters = [
+    ["radius",         "Ang",       20, [0, inf],    "volume",      "Cylinder core radius"],
+    ["rim_thickness",  "Ang",       10, [0, inf],    "volume",      "Rim shell thickness"],
+    ["face_thickness", "Ang",       10, [0, inf],    "volume",      "Cylinder face thickness"],
+    ["length",         "Ang",      400, [0, inf],    "volume",      "Cylinder length"],
+    ["core_sld",       "1e-6/Ang^2", 1, [-inf, inf], "",            "Cylinder core scattering length density"],
+    ["face_sld",       "1e-6/Ang^2", 4, [-inf, inf], "",            "Cylinder face scattering length density"],
+    ["rim_sld",        "1e-6/Ang^2", 4, [-inf, inf], "",            "Cylinder rim scattering length density"],
+    ["solvent_sld",    "1e-6/Ang^2", 1, [-inf, inf], "",            "Solvent scattering length density"],
+    ["theta",          "degrees",   90, [-inf, inf], "orientation", "In plane angle"],
+    ["phi",            "degrees",    0, [-inf, inf], "orientation", "Out of plane angle"],
+    ]
+
+# pylint: enable=bad-whitespace, line-too-long
+
+source = ["lib/J1.c", "lib/gauss76.c", "core_shell_bicelle.c"]
+
+demo = dict(scale=1, background=0,
+            radius=20.0,
+            rim_thickness=10.0,
+            face_thickness=10.0,
+            length=400.0,
+            core_sld=1.0,
+            face_sld=4.0,
+            rim_sld=4.0,
+            solvent_sld=1.0,
+            theta=90,
+            phi=0)
+
+oldname = 'CoreShellBicelleModel'
+
+oldpars = dict(rim_thickness='rim_thick',
+               face_thickness='face_thick',
+               theta='axis_theta',
+               phi='axis_phi')
+
+qx, qy = 0.4 * cos(90), 0.5 * sin(0)
+tests = [
+    # Accuracy tests based on content in test/utest_other_models.py
+    [{'radius': 20.0,
+      'rim_thickness': 10.0,
+      'face_thickness': 10.0,
+      'length': 400.0,
+      'core_sld': 1.0,
+      'face_sld': 4.0,
+      'rim_sld': 4.0,
+      'solvent_sld': 1.0,
+      'background': 0.0,
+     }, 0.001, 353.550],
+
+    [{'radius': 20.0,
+      'rim_thickness': 10.0,
+      'face_thickness': 10.0,
+      'length': 400.0,
+      'core_sld': 1.0,
+      'face_sld': 4.0,
+      'rim_sld': 4.0,
+      'solvent_sld': 1.0,
+      'theta': 90.0,
+      'phi': 0.0,
+      'background': 0.00,
+     }, (qx, qy), 24.9167],
+
+    # Additional tests with larger range of parameters
+    [{'radius': 3.0,
+      'rim_thickness': 100.0,
+      'face_thickness': 100.0,
+      'length': 1200.0,
+      'core_sld': 5.0,
+      'face_sld': 41.0,
+      'rim_sld': 42.0,
+      'solvent_sld': 21.0,
+     }, 0.05, 1670.1828],
+    ]
