diff --git a/GPy/kern/__init__.py b/GPy/kern/__init__.py
index 97c1d88f..6d9cf07a 100644
--- a/GPy/kern/__init__.py
+++ b/GPy/kern/__init__.py
@@ -1,8 +1,7 @@
 # Copyright (c) 2012, GPy authors (see AUTHORS.txt).
 # Licensed under the BSD 3-clause license (see LICENSE.txt)
 
-
-from constructors import rbf, Matern32, Matern52, exponential, linear, white, bias, finite_dimensional, spline, Brownian, periodic_exponential, periodic_Matern32, periodic_Matern52, prod, symmetric, Coregionalise, rational_quadratic, Fixed, rbfcos, IndependentOutputs
+from constructors import *
 try:
     from constructors import rbf_sympy, sympykern # these depend on sympy
 except:
diff --git a/GPy/kern/constructors.py b/GPy/kern/constructors.py
index e2c21f15..697f3554 100644
--- a/GPy/kern/constructors.py
+++ b/GPy/kern/constructors.py
@@ -1,33 +1,9 @@
 # Copyright (c) 2012, GPy authors (see AUTHORS.txt).
 # Licensed under the BSD 3-clause license (see LICENSE.txt)
 
-
 import numpy as np
 from kern import kern
-
-from rbf import rbf as rbfpart
-from white import white as whitepart
-from linear import linear as linearpart
-from exponential import exponential as exponentialpart
-from Matern32 import Matern32 as Matern32part
-from Matern52 import Matern52 as Matern52part
-from bias import bias as biaspart
-from fixed import Fixed as fixedpart
-from finite_dimensional import finite_dimensional as finite_dimensionalpart
-from spline import spline as splinepart
-from Brownian import Brownian as Brownianpart
-from periodic_exponential import periodic_exponential as periodic_exponentialpart
-from periodic_Matern32 import periodic_Matern32 as periodic_Matern32part
-from periodic_Matern52 import periodic_Matern52 as periodic_Matern52part
-from prod import prod as prodpart
-from symmetric import symmetric as symmetric_part
-from coregionalise import Coregionalise as coregionalise_part
-from rational_quadratic import rational_quadratic as rational_quadraticpart
-from rbfcos import rbfcos as rbfcospart
-from independent_outputs import IndependentOutputs as independent_output_part
-#TODO these s=constructors are not as clean as we'd like. Tidy the code up
-#using meta-classes to make the objects construct properly wthout them.
-
+import parts
 
 def rbf(input_dim,variance=1., lengthscale=None,ARD=False):
     """
@@ -42,7 +18,7 @@ def rbf(input_dim,variance=1., lengthscale=None,ARD=False):
     :param ARD: Auto Relevance Determination (one lengthscale per dimension)
     :type ARD: Boolean
     """
-    part = rbfpart(input_dim,variance,lengthscale,ARD)
+    part = parts.rbf.RBF(input_dim,variance,lengthscale,ARD)
     return kern(input_dim, [part])
 
 def linear(input_dim,variances=None,ARD=False):
@@ -55,7 +31,7 @@ def linear(input_dim,variances=None,ARD=False):
      variances (np.ndarray)
      ARD (boolean)
     """
-    part = linearpart(input_dim,variances,ARD)
+    part = parts.linear.Linear(input_dim,variances,ARD)
     return kern(input_dim, [part])
 
 def white(input_dim,variance=1.):
@@ -67,7 +43,7 @@ def white(input_dim,variance=1.):
     input_dimD (int), obligatory
      variance (float)
     """
-    part = whitepart(input_dim,variance)
+    part = parts.white.White(input_dim,variance)
     return kern(input_dim, [part])
 
 def exponential(input_dim,variance=1., lengthscale=None, ARD=False):
@@ -83,7 +59,7 @@ def exponential(input_dim,variance=1., lengthscale=None, ARD=False):
     :param ARD: Auto Relevance Determination (one lengthscale per dimension)
     :type ARD: Boolean
     """
-    part = exponentialpart(input_dim,variance, lengthscale, ARD)
+    part = parts.exponential.Exponential(input_dim,variance, lengthscale, ARD)
     return kern(input_dim, [part])
 
 def Matern32(input_dim,variance=1., lengthscale=None, ARD=False):
@@ -99,7 +75,7 @@ def Matern32(input_dim,variance=1., lengthscale=None, ARD=False):
     :param ARD: Auto Relevance Determination (one lengthscale per dimension)
     :type ARD: Boolean
     """
-    part = Matern32part(input_dim,variance, lengthscale, ARD)
+    part = parts.Matern32.Matern32(input_dim,variance, lengthscale, ARD)
     return kern(input_dim, [part])
 
 def Matern52(input_dim, variance=1., lengthscale=None, ARD=False):
@@ -115,7 +91,7 @@ def Matern52(input_dim, variance=1., lengthscale=None, ARD=False):
     :param ARD: Auto Relevance Determination (one lengthscale per dimension)
     :type ARD: Boolean
     """
-    part = Matern52part(input_dim, variance, lengthscale, ARD)
+    part = parts.Matern52.Matern52(input_dim, variance, lengthscale, ARD)
     return kern(input_dim, [part])
 
 def bias(input_dim, variance=1.):
@@ -127,7 +103,7 @@ def bias(input_dim, variance=1.):
      input_dim (int), obligatory
      variance (float)
     """
-    part = biaspart(input_dim, variance)
+    part = parts.bias.Bias(input_dim, variance)
     return kern(input_dim, [part])
 
 def finite_dimensional(input_dim, F, G, variances=1., weights=None):
@@ -138,7 +114,7 @@ def finite_dimensional(input_dim, F, G, variances=1., weights=None):
     G: np.array with shape (n,n) - the Gram matrix associated to F
     variances : np.ndarray with shape (n,)
     """
-    part = finite_dimensionalpart(input_dim, F, G, variances, weights)
+    part = parts.finite_dimensional.FiniteDimensional(input_dim, F, G, variances, weights)
     return kern(input_dim, [part])
 
 def spline(input_dim, variance=1.):
@@ -150,7 +126,7 @@ def spline(input_dim, variance=1.):
     :param variance: the variance of the kernel
     :type variance: float
     """
-    part = splinepart(input_dim, variance)
+    part = parts.spline.Spline(input_dim, variance)
     return kern(input_dim, [part])
 
 def Brownian(input_dim, variance=1.):
@@ -162,7 +138,7 @@ def Brownian(input_dim, variance=1.):
     :param variance: the variance of the kernel
     :type variance: float
     """
-    part = Brownianpart(input_dim, variance)
+    part = parts.Brownian.Brownian(input_dim, variance)
     return kern(input_dim, [part])
 
 try:
@@ -215,7 +191,7 @@ def periodic_exponential(input_dim=1, variance=1., lengthscale=None, period=2 *
     :param n_freq: the number of frequencies considered for the periodic subspace
     :type n_freq: int
     """
-    part = periodic_exponentialpart(input_dim, variance, lengthscale, period, n_freq, lower, upper)
+    part = parts.periodic_exponential.PeriodicExponential(input_dim, variance, lengthscale, period, n_freq, lower, upper)
     return kern(input_dim, [part])
 
 def periodic_Matern32(input_dim, variance=1., lengthscale=None, period=2 * np.pi, n_freq=10, lower=0., upper=4 * np.pi):
@@ -233,7 +209,7 @@ def periodic_Matern32(input_dim, variance=1., lengthscale=None, period=2 * np.pi
      :param n_freq: the number of frequencies considered for the periodic subspace
      :type n_freq: int
     """
-    part = periodic_Matern32part(input_dim, variance, lengthscale, period, n_freq, lower, upper)
+    part = parts.periodic_Matern32.PeriodicMatern32(input_dim, variance, lengthscale, period, n_freq, lower, upper)
     return kern(input_dim, [part])
 
 def periodic_Matern52(input_dim, variance=1., lengthscale=None, period=2 * np.pi, n_freq=10, lower=0., upper=4 * np.pi):
@@ -251,7 +227,7 @@ def periodic_Matern52(input_dim, variance=1., lengthscale=None, period=2 * np.pi
      :param n_freq: the number of frequencies considered for the periodic subspace
      :type n_freq: int
     """
-    part = periodic_Matern52part(input_dim, variance, lengthscale, period, n_freq, lower, upper)
+    part = parts.periodic_Matern52part(input_dim, variance, lengthscale, period, n_freq, lower, upper)
     return kern(input_dim, [part])
 
 def prod(k1,k2,tensor=False):
@@ -262,7 +238,7 @@ def prod(k1,k2,tensor=False):
     :type k1, k2: kernpart
     :rtype: kernel object
     """
-    part = prodpart(k1,k2,tensor)
+    part = parts.prodpart(k1,k2,tensor)
     return kern(part.input_dim, [part])
 
 def symmetric(k):
@@ -270,11 +246,11 @@ def symmetric(k):
     Construct a symmetrical kernel from an existing kernel
     """
     k_ = k.copy()
-    k_.parts = [symmetric_part(p) for p in k.parts]
+    k_.parts = [symmetric.Symmetric(p) for p in k.parts]
     return k_
 
-def Coregionalise(Nout,R=1, W=None, kappa=None):
-    p = coregionalise_part(Nout,R,W,kappa)
+def coregionalise(Nout,R=1, W=None, kappa=None):
+    p = parts.coregionalise.Coregionalise(Nout,R,W,kappa)
     return kern(1,[p])
 
 
@@ -291,10 +267,10 @@ def rational_quadratic(input_dim, variance=1., lengthscale=1., power=1.):
     :rtype: kern object
 
     """
-    part = rational_quadraticpart(input_dim, variance, lengthscale, power)
+    part = parts.rational_quadratic.RationalQuadratic(input_dim, variance, lengthscale, power)
     return kern(input_dim, [part])
 
-def Fixed(input_dim, K, variance=1.):
+def fixed(input_dim, K, variance=1.):
     """
      Construct a Fixed effect kernel.
 
@@ -304,23 +280,21 @@ def Fixed(input_dim, K, variance=1.):
      K (np.array), obligatory
      variance (float)
     """
-    part = fixedpart(input_dim, K, variance)
+    part = parts.fixed.Fixed(input_dim, K, variance)
     return kern(input_dim, [part])
 
 def rbfcos(input_dim, variance=1., frequencies=None, bandwidths=None, ARD=False):
     """
     construct a rbfcos kernel
     """
-    part = rbfcospart(input_dim, variance, frequencies, bandwidths, ARD)
+    part = parts.rbfcos.RBFCos(input_dim, variance, frequencies, bandwidths, ARD)
     return kern(input_dim, [part])
 
-def IndependentOutputs(k):
+def independent_outputs(k):
     """
     Construct a kernel with independent outputs from an existing kernel
     """
     for sl in k.input_slices:
         assert (sl.start is None) and (sl.stop is None), "cannot adjust input slices! (TODO)"
-    parts = [independent_output_part(p) for p in k.parts]
+    parts = [independent_outputs.IndependentOutputs(p) for p in k.parts]
     return kern(k.input_dim+1,parts)
-
-
diff --git a/GPy/kern/kern.py b/GPy/kern/kern.py
index 76a0d99e..ad3ed3f9 100644
--- a/GPy/kern/kern.py
+++ b/GPy/kern/kern.py
@@ -1,13 +1,12 @@
 # Copyright (c) 2012, GPy authors (see AUTHORS.txt).
 # Licensed under the BSD 3-clause license (see LICENSE.txt)
 
-
 import numpy as np
 import pylab as pb
 from ..core.parameterised import Parameterised
-from kernpart import Kernpart
+from parts.kernpart import Kernpart
 import itertools
-from prod import prod
+from parts.prod import Prod as prod
 
 class kern(Parameterised):
     def __init__(self, input_dim, parts=[], input_slices=None):
diff --git a/GPy/kern/Brownian.py b/GPy/kern/parts/Brownian.py
similarity index 100%
rename from GPy/kern/Brownian.py
rename to GPy/kern/parts/Brownian.py
diff --git a/GPy/kern/Matern32.py b/GPy/kern/parts/Matern32.py
similarity index 100%
rename from GPy/kern/Matern32.py
rename to GPy/kern/parts/Matern32.py
diff --git a/GPy/kern/Matern52.py b/GPy/kern/parts/Matern52.py
similarity index 100%
rename from GPy/kern/Matern52.py
rename to GPy/kern/parts/Matern52.py
diff --git a/GPy/kern/bias.py b/GPy/kern/parts/bias.py
similarity index 99%
rename from GPy/kern/bias.py
rename to GPy/kern/parts/bias.py
index 8ec3741d..2b72e7c9 100644
--- a/GPy/kern/bias.py
+++ b/GPy/kern/parts/bias.py
@@ -6,7 +6,7 @@ from kernpart import Kernpart
 import numpy as np
 import hashlib
 
-class bias(Kernpart):
+class Bias(Kernpart):
     def __init__(self,input_dim,variance=1.):
         """
         :param input_dim: the number of input dimensions
diff --git a/GPy/kern/coregionalise.py b/GPy/kern/parts/coregionalise.py
similarity index 100%
rename from GPy/kern/coregionalise.py
rename to GPy/kern/parts/coregionalise.py
diff --git a/GPy/kern/exponential.py b/GPy/kern/parts/exponential.py
similarity index 99%
rename from GPy/kern/exponential.py
rename to GPy/kern/parts/exponential.py
index 5cb0f584..d8cf76f7 100644
--- a/GPy/kern/exponential.py
+++ b/GPy/kern/parts/exponential.py
@@ -6,7 +6,7 @@ from kernpart import Kernpart
 import numpy as np
 from scipy import integrate
 
-class exponential(Kernpart):
+class Exponential(Kernpart):
     """
     Exponential kernel (aka Ornstein-Uhlenbeck or Matern 1/2)
 
diff --git a/GPy/kern/finite_dimensional.py b/GPy/kern/parts/finite_dimensional.py
similarity index 97%
rename from GPy/kern/finite_dimensional.py
rename to GPy/kern/parts/finite_dimensional.py
index b23ddb16..6cc2325f 100644
--- a/GPy/kern/finite_dimensional.py
+++ b/GPy/kern/parts/finite_dimensional.py
@@ -4,9 +4,9 @@
 
 from kernpart import Kernpart
 import numpy as np
-from ..util.linalg import pdinv,mdot
+from ...util.linalg import pdinv,mdot
 
-class finite_dimensional(Kernpart):
+class FiniteDimensional(Kernpart):
     def __init__(self, input_dim, F, G, variance=1., weights=None):
         """
         Argumnents
diff --git a/GPy/kern/fixed.py b/GPy/kern/parts/fixed.py
similarity index 100%
rename from GPy/kern/fixed.py
rename to GPy/kern/parts/fixed.py
diff --git a/GPy/kern/independent_outputs.py b/GPy/kern/parts/independent_outputs.py
similarity index 100%
rename from GPy/kern/independent_outputs.py
rename to GPy/kern/parts/independent_outputs.py
diff --git a/GPy/kern/kernpart.py b/GPy/kern/parts/kernpart.py
similarity index 100%
rename from GPy/kern/kernpart.py
rename to GPy/kern/parts/kernpart.py
diff --git a/GPy/kern/linear.py b/GPy/kern/parts/linear.py
similarity index 99%
rename from GPy/kern/linear.py
rename to GPy/kern/parts/linear.py
index 0ed887d2..c1da3944 100644
--- a/GPy/kern/linear.py
+++ b/GPy/kern/parts/linear.py
@@ -4,10 +4,10 @@
 
 from kernpart import Kernpart
 import numpy as np
-from ..util.linalg import tdot
+from ...util.linalg import tdot
 from scipy import weave
 
-class linear(Kernpart):
+class Linear(Kernpart):
     """
     Linear kernel
 
diff --git a/GPy/kern/periodic_Matern32.py b/GPy/kern/parts/periodic_Matern32.py
similarity index 99%
rename from GPy/kern/periodic_Matern32.py
rename to GPy/kern/parts/periodic_Matern32.py
index 664a5183..5693085d 100644
--- a/GPy/kern/periodic_Matern32.py
+++ b/GPy/kern/parts/periodic_Matern32.py
@@ -7,7 +7,7 @@ import numpy as np
 from GPy.util.linalg import mdot
 from GPy.util.decorators import silence_errors
 
-class periodic_Matern32(Kernpart):
+class PeriodicMatern32(Kernpart):
     """
     Kernel of the periodic subspace (up to a given frequency) of a Matern 3/2 RKHS. Only defined for input_dim=1.
 
diff --git a/GPy/kern/periodic_Matern52.py b/GPy/kern/parts/periodic_Matern52.py
similarity index 99%
rename from GPy/kern/periodic_Matern52.py
rename to GPy/kern/parts/periodic_Matern52.py
index c01d2d26..7b5ae846 100644
--- a/GPy/kern/periodic_Matern52.py
+++ b/GPy/kern/parts/periodic_Matern52.py
@@ -7,7 +7,7 @@ import numpy as np
 from GPy.util.linalg import mdot
 from GPy.util.decorators import silence_errors
 
-class periodic_Matern52(Kernpart):
+class PeriodicMatern52(Kernpart):
     """
     Kernel of the periodic subspace (up to a given frequency) of a Matern 5/2 RKHS. Only defined for input_dim=1.
 
diff --git a/GPy/kern/periodic_exponential.py b/GPy/kern/parts/periodic_exponential.py
similarity index 99%
rename from GPy/kern/periodic_exponential.py
rename to GPy/kern/parts/periodic_exponential.py
index fcaf6420..36b7b9ac 100644
--- a/GPy/kern/periodic_exponential.py
+++ b/GPy/kern/parts/periodic_exponential.py
@@ -7,7 +7,7 @@ import numpy as np
 from GPy.util.linalg import mdot
 from GPy.util.decorators import silence_errors
 
-class periodic_exponential(Kernpart):
+class PeriodicExponential(Kernpart):
     """
     Kernel of the periodic subspace (up to a given frequency) of a exponential (Matern 1/2) RKHS. Only defined for input_dim=1.
 
diff --git a/GPy/kern/prod.py b/GPy/kern/parts/prod.py
similarity index 99%
rename from GPy/kern/prod.py
rename to GPy/kern/parts/prod.py
index 493129c6..db31c626 100644
--- a/GPy/kern/prod.py
+++ b/GPy/kern/parts/prod.py
@@ -5,7 +5,7 @@ from kernpart import Kernpart
 import numpy as np
 import hashlib
 
-class prod(Kernpart):
+class Prod(Kernpart):
     """
     Computes the product of 2 kernels
 
diff --git a/GPy/kern/prod_orthogonal.py b/GPy/kern/parts/prod_orthogonal.py
similarity index 100%
rename from GPy/kern/prod_orthogonal.py
rename to GPy/kern/parts/prod_orthogonal.py
diff --git a/GPy/kern/rational_quadratic.py b/GPy/kern/parts/rational_quadratic.py
similarity index 98%
rename from GPy/kern/rational_quadratic.py
rename to GPy/kern/parts/rational_quadratic.py
index d1e7a7e3..d92b43db 100644
--- a/GPy/kern/rational_quadratic.py
+++ b/GPy/kern/parts/rational_quadratic.py
@@ -5,7 +5,7 @@
 from kernpart import Kernpart
 import numpy as np
 
-class rational_quadratic(Kernpart):
+class RationalQuadratic(Kernpart):
     """
     rational quadratic kernel
 
diff --git a/GPy/kern/rbf.py b/GPy/kern/parts/rbf.py
similarity index 99%
rename from GPy/kern/rbf.py
rename to GPy/kern/parts/rbf.py
index 5686d7bd..aa07e15d 100644
--- a/GPy/kern/rbf.py
+++ b/GPy/kern/parts/rbf.py
@@ -6,9 +6,9 @@ from kernpart import Kernpart
 import numpy as np
 import hashlib
 from scipy import weave
-from ..util.linalg import tdot
+from ...util.linalg import tdot
 
-class rbf(Kernpart):
+class RBF(Kernpart):
     """
     Radial Basis Function kernel, aka squared-exponential, exponentiated quadratic or Gaussian kernel:
 
diff --git a/GPy/kern/rbfcos.py b/GPy/kern/parts/rbfcos.py
similarity index 99%
rename from GPy/kern/rbfcos.py
rename to GPy/kern/parts/rbfcos.py
index b1e99d3c..4d09dfdb 100644
--- a/GPy/kern/rbfcos.py
+++ b/GPy/kern/parts/rbfcos.py
@@ -6,7 +6,7 @@
 from kernpart import Kernpart
 import numpy as np
 
-class rbfcos(Kernpart):
+class RBFCos(Kernpart):
     def __init__(self,input_dim,variance=1.,frequencies=None,bandwidths=None,ARD=False):
         self.input_dim = input_dim
         self.name = 'rbfcos'
diff --git a/GPy/kern/spline.py b/GPy/kern/parts/spline.py
similarity index 98%
rename from GPy/kern/spline.py
rename to GPy/kern/parts/spline.py
index f2802180..e675420c 100644
--- a/GPy/kern/spline.py
+++ b/GPy/kern/parts/spline.py
@@ -9,7 +9,7 @@ def theta(x):
     """Heaviside step function"""
     return np.where(x>=0.,1.,0.)
 
-class spline(Kernpart):
+class Spline(Kernpart):
     """
     Spline kernel
 
diff --git a/GPy/kern/symmetric.py b/GPy/kern/parts/symmetric.py
similarity index 99%
rename from GPy/kern/symmetric.py
rename to GPy/kern/parts/symmetric.py
index c7099a6f..bbdd5ac0 100644
--- a/GPy/kern/symmetric.py
+++ b/GPy/kern/parts/symmetric.py
@@ -4,7 +4,7 @@
 from kernpart import Kernpart
 import numpy as np
 
-class symmetric(Kernpart):
+class Symmetric(Kernpart):
     """
     Symmetrical kernels
 
diff --git a/GPy/kern/sympy_helpers.cpp b/GPy/kern/parts/sympy_helpers.cpp
similarity index 100%
rename from GPy/kern/sympy_helpers.cpp
rename to GPy/kern/parts/sympy_helpers.cpp
diff --git a/GPy/kern/sympy_helpers.h b/GPy/kern/parts/sympy_helpers.h
similarity index 100%
rename from GPy/kern/sympy_helpers.h
rename to GPy/kern/parts/sympy_helpers.h
diff --git a/GPy/kern/sympykern.py b/GPy/kern/parts/sympykern.py
similarity index 100%
rename from GPy/kern/sympykern.py
rename to GPy/kern/parts/sympykern.py
diff --git a/GPy/kern/white.py b/GPy/kern/parts/white.py
similarity index 98%
rename from GPy/kern/white.py
rename to GPy/kern/parts/white.py
index 41f075c3..8241bf55 100644
--- a/GPy/kern/white.py
+++ b/GPy/kern/parts/white.py
@@ -1,10 +1,10 @@
 # Copyright (c) 2012, GPy authors (see AUTHORS.txt).
 # Licensed under the BSD 3-clause license (see LICENSE.txt)
 
-
 from kernpart import Kernpart
 import numpy as np
-class white(Kernpart):
+
+class White(Kernpart):
     """
     White noise kernel.
 
diff --git a/GPy/testing/kernel_tests.py b/GPy/testing/kernel_tests.py
index 98c75827..8a4f8b16 100644
--- a/GPy/testing/kernel_tests.py
+++ b/GPy/testing/kernel_tests.py
@@ -21,7 +21,7 @@ class KernelTests(unittest.TestCase):
         """
         X = np.random.rand(30, 4)
         K = np.dot(X, X.T)
-        kernel = GPy.kern.Fixed(4, K)
+        kernel = GPy.kern.fixed(4, K)
         Y = np.ones((30,1))
         m = GPy.models.GPRegression(X,Y,kernel=kernel)
         self.assertTrue(m.checkgrad())
@@ -36,7 +36,7 @@ class KernelTests(unittest.TestCase):
         Y = np.vstack((Y1,Y2))
 
         k1 = GPy.kern.rbf(1) + GPy.kern.bias(1)
-        k2 = GPy.kern.Coregionalise(2,1)
+        k2 = GPy.kern.coregionalise(2,1)
         k = k1.prod(k2,tensor=True)
         m = GPy.models.GPRegression(X,Y,kernel=k)
         self.assertTrue(m.checkgrad())