diff --git a/src/maths/types.nim b/src/maths/types.nim
index 7fc84c2..42c56d6 100644
--- a/src/maths/types.nim
+++ b/src/maths/types.nim
@@ -267,7 +267,7 @@ template declaredVector*(x:ToSingle):untyped = isVector(x[])
 template declaredMatrix*(x:ToSingle):untyped = isMatrix(x[])
 template re*(x: ToSingle): untyped = toSingle(x[].re)
 template im*(x: ToSingle): untyped = toSingle(x[].im)
-template simdType*(x: ToSingle): untyped = simdType(x[])
+#template simdType*(x: ToSingle): untyped = simdType(x[])
 template numberType*(x: ToSingle): untyped = float32
 
 type
@@ -328,7 +328,7 @@ template declaredVector*(x:ToDouble):untyped = isVector(x[])
 template declaredMatrix*(x:ToDouble):untyped = isMatrix(x[])
 template re*(x:ToDouble):untyped = toDouble(x[].re)
 template im*(x:ToDouble):untyped = toDouble(x[].im)
-template simdType*(x: ToDouble): untyped = simdType(x[])
+#template simdType*(x: ToDouble): untyped = simdType(x[])
 #macro dump2(x: typed): auto =
 #  result = newEmptyNode()
 #  echo x.treerepr
@@ -505,7 +505,7 @@ forwardFunc(Indexed, norm2)
 template len*(x:Indexed):untyped = obj(x).len  # FIXME Simd types can use Indexed
 template nrows*(x:Indexed):untyped = obj(x).nrows
 template ncols*(x:Indexed):untyped = obj(x).ncols
-template simdType*(x: Indexed): untyped = simdType(obj(x))
+#template simdType*(x: Indexed): untyped = simdType(obj(x))
 template numberType*(x: Indexed): untyped = numberType(obj(x))
 #template declaredVector*(x:Indexed):untyped = isVector(x.obj)
 #template declaredMatrix*(x:Indexed):untyped = isMatrix(x.obj)
diff --git a/src/physics/colorTensor.nim b/src/physics/colorTensor.nim
index e9d8bfa..306baf8 100644
--- a/src/physics/colorTensor.nim
+++ b/src/physics/colorTensor.nim
@@ -8,6 +8,7 @@ type
   Color3*[T] = Color[T]
   Color4*[T] = Color[T]
 
+proc has*(t: typedesc[ColorType], x: typedesc): bool = x is ColorType
 #[
 template getNc*[K,T](x: SomeTensor[K,T]): int =
   when K is ColorType:
diff --git a/src/physics/spin.nim b/src/physics/spin.nim
index 2a07f1b..fc8ef84 100644
--- a/src/physics/spin.nim
+++ b/src/physics/spin.nim
@@ -1,5 +1,5 @@
-import spinOld
-export spinOld
+#import spinOld
+#export spinOld
 
-#import spinTensor
-#export spinTensor
+import spinTensor
+export spinTensor
diff --git a/src/physics/spinTensor.nim b/src/physics/spinTensor.nim
new file mode 100644
index 0000000..2841256
--- /dev/null
+++ b/src/physics/spinTensor.nim
@@ -0,0 +1,306 @@
+import tensorwrap
+export tensorwrap
+import maths/[complexType, matrixConcept]
+
+type
+  SpinType* = object
+  Spin*[T] = TensorObj[SpinType,T]
+  Spin2*[T] = Spin[T]
+  Spin3*[T] = Spin[T]
+  Spin4*[T] = Spin[T]
+template asSpin*[T](x: typedesc[T]): typedesc = Spin[T]
+template asSpin*[T](x: T): auto = Spin[T](obj: x)
+
+template spinVector*(x: static[int], a: untyped): untyped =
+  #const
+  #  I:int = x
+  #type
+  #  E = T
+  #  VA = VectorArray[I,E]
+  #  VAO = VectorArrayObj[I,E]
+  #Spin[VA](v: VA(v: a))
+  #Spin[VA](v: VA(v: VAO(vec: a)))
+  #asSpin(VA(v: a))
+  #static: echo "spinVector"
+  asSpin(asVectorArray(a))
+  #let t1 = asVectorArray(a)
+  #static: echo "spinVector1"
+  #let t = asSpin(t1)
+  #static: echo "spinVector2"
+  #t
+template spinVector*[T](x: static[int], a: typedesc[T]): untyped =
+  asSpin(asVectorArray(x, type T))
+template spinMatrix*[T](x,y:static[int], a: untyped): untyped =
+  const
+    I:int = x
+    J:int = y
+  type
+    E = T
+    #MA = MatrixArray[I,J,E]
+    MAO = MatrixArrayObj[I,J,E]
+  asSpin(asMatrix(MAO(mat: a)))
+#template spinMatrix*[I,J:static[int],T](a: untyped): untyped =
+#  Spin[MatrixArray[I,J,T]](v: MatrixArray[I,J,T](v: MatrixArrayObj[I,J,T](mat: a)))
+#template spinMatrix*(I,J,T,a: untyped): untyped =
+#  Spin[MatrixArray[I,J,T]](v: MatrixArray[I,J,T](v: MatrixArrayObj[I,J,T](mat: a)))
+
+const z0 = newComplex( 0.0,  0.0)
+const z1 = newComplex( 1.0,  0.0)
+const zi = newComplex( 0.0,  1.0)
+const n1 = newComplex(-1.0,  0.0)
+const ni = newComplex( 0.0, -1.0)
+
+template s(r,c,x: untyped): untyped =
+  spinMatrix[ComplexType[float]](r,c,x)
+template g(x: untyped): untyped = s(4,4,x)
+template p(x: untyped): untyped = s(2,4,x)
+template r(x: untyped): untyped = s(4,2,x)
+
+const
+  gamma0* = g([[ z1, z0, z0, z0 ],
+               [ z0, z1, z0, z0 ],
+               [ z0, z0, z1, z0 ],
+               [ z0, z0, z0, z1 ]])
+  gamma1* = g([[ z0, z0, z0, zi ],
+               [ z0, z0, zi, z0 ],
+               [ z0, ni, z0, z0 ],
+               [ ni, z0, z0, z0 ]])
+  gamma2* = g([[ z0, z0, z0, n1 ],
+               [ z0, z0, z1, z0 ],
+               [ z0, z1, z0, z0 ],
+               [ n1, z0, z0, z0 ]])
+  gamma3* = g([[ z0, z0, zi, z0 ],
+               [ z0, z0, z0, ni ],
+               [ ni, z0, z0, z0 ],
+               [ z0, zi, z0, z0 ]])
+  gamma4* = g([[ z0, z0, z1, z0 ],
+               [ z0, z0, z0, z1 ],
+               [ z1, z0, z0, z0 ],
+               [ z0, z1, z0, z0 ]])
+  gamma5* = g([[ z1, z0, z0, z0 ],
+               [ z0, z1, z0, z0 ],
+               [ z0, z0, n1, z0 ],
+               [ z0, z0, z0, n1 ]])
+
+  spprojmat1p* = p([[ z1, z0, z0, zi ],
+                    [ z0, z1, zi, z0 ]])
+  spprojmat1m* = p([[ z1, z0, z0, ni ],
+                    [ z0, z1, ni, z0 ]])
+  spprojmat2p* = p([[ z1, z0, z0, n1 ],
+                    [ z0, z1, z1, z0 ]])
+  spprojmat2m* = p([[ z1, z0, z0, z1 ],
+                    [ z0, z1, n1, z0 ]])
+  spprojmat3p* = p([[ z1, z0, zi, z0 ],
+                    [ z0, z1, z0, ni ]])
+  spprojmat3m* = p([[ z1, z0, ni, z0 ],
+                    [ z0, z1, z0, zi ]])
+  spprojmat4p* = p([[ z1, z0, z1, z0 ],
+                    [ z0, z1, z0, z1 ]])
+  spprojmat4m* = p([[ z1, z0, n1, z0 ],
+                    [ z0, z1, z0, n1 ]])
+
+  spreconmat1p* = r([[ z1, z0 ],
+                     [ z0, z1 ],
+                     [ z0, ni ],
+                     [ ni, z0 ]])
+  spreconmat1m* = r([[ z1, z0 ],
+                     [ z0, z1 ],
+                     [ z0, zi ],
+                     [ zi, z0 ]])
+  spreconmat2p* = r([[ z1, z0 ],
+                     [ z0, z1 ],
+                     [ z0, z1 ],
+                     [ n1, z0 ]])
+  spreconmat2m* = r([[ z1, z0 ],
+                     [ z0, z1 ],
+                     [ z0, n1 ],
+                     [ z1, z0 ]])
+  spreconmat3p* = r([[ z1, z0 ],
+                     [ z0, z1 ],
+                     [ ni, z0 ],
+                     [ z0, zi ]])
+  spreconmat3m* = r([[ z1, z0 ],
+                     [ z0, z1 ],
+                     [ zi, z0 ],
+                     [ z0, ni ]])
+  spreconmat4p* = r([[ z1, z0 ],
+                     [ z0, z1 ],
+                     [ z1, z0 ],
+                     [ z0, z1 ]])
+  spreconmat4m* = r([[ z1, z0 ],
+                     [ z0, z1 ],
+                     [ n1, z0 ],
+                     [ z0, n1 ]])
+
+template I(x: typed): untyped =
+  newImag(1)*x
+template mI(x: typed): untyped =
+  newImag(-1)*x
+
+proc spproj1p*(r: var auto, x: auto) =
+  ## r: HalfFermion
+  ## x: DiracFermion
+  #let nc = x[0].len
+  #for i in 0..<nc:
+  #  r[0][i] = x[0][i] + x[2][i]
+  #  r[1][i] = x[1][i] + x[3][i]
+  r := spprojmat1p * x
+
+#[
+template spproj1p*(x: auto): untyped = spprojmat1p * x
+template spproj2p*(x: auto): untyped = spprojmat2p * x
+template spproj3p*(x: auto): untyped = spprojmat3p * x
+template spproj4p*(x: auto): untyped = spprojmat4p * x
+template spproj1m*(x: auto): untyped = spprojmat1m * x
+template spproj2m*(x: auto): untyped = spprojmat2m * x
+template spproj3m*(x: auto): untyped = spprojmat3m * x
+template spproj4m*(x: auto): untyped = spprojmat4m * x
+]#
+
+#template spproj1pU*(x: Spin): untyped =
+  #let v0 = x[][0] + I(x[][3])
+  #let v1 = x[][1] + I(x[][2])
+  #spinVector[type(v0)](2,[v0,v1])
+  #spinVector[type(x[][0])](2,[x[][0]+I(x[][3]),x[][1]+I(x[][2])])
+#  spinVector(2,[x[][0]+I(x[][3]),x[][1]+I(x[][2])])
+#template spproj1p*(x: Spin): untyped =
+#  flattenCallArgs(spproj1pU, x)
+template spproj1p*(xx: Spin): untyped =
+  let xp = getPtr xx; template x:untyped {.gensym.} = xp[]
+  #let v = [ x[0]+I(x[3]), x[1]+I(x[2]) ]
+  #spinVector(2, v)
+  var v {.noInit.}: spinVector(2,evalType(x[0]))
+  #static: echo "spproj1p: ", v.type
+  #v[0] = x[0]+I(x[3])
+  #v[1] = x[1]+I(x[2])
+  let t0 = I(x[3])
+  add(v[0], x[0], t0)
+  let t1 = I(x[2])
+  add(v[1], x[1], t1)
+  v
+template spproj2p*(xx: Spin): untyped =
+  let xp = getPtr xx; template x:untyped {.gensym.} = xp[]
+  let v = [ x[0]-x[3], x[1]+x[2] ]
+  spinVector(2, v)
+template spproj3p*(xx: Spin): untyped =
+  let xp = getPtr xx; template x:untyped {.gensym.} = xp[]
+  let v = [ x[0]+I(x[2]), x[1]-I(x[3]) ]
+  spinVector(2, v)
+template spproj4p*(xx: Spin): untyped =
+  let xp = getPtr xx; template x:untyped {.gensym.} = xp[]
+  let v = [ x[0]+x[2], x[1]+x[3] ]
+  spinVector(2, v)
+template spproj1m*(xx: Spin): untyped =
+  let xp = getPtr xx; template x:untyped {.gensym.} = xp[]
+  let v = [ x[0]-I(x[3]), x[1]-I(x[2]) ]
+  spinVector(2, v)
+template spproj2m*(xx: Spin): untyped =
+  let xp = getPtr xx; template x:untyped {.gensym.} = xp[]
+  let v = [ x[0]+x[3], x[1]-x[2] ]
+  spinVector(2, v)
+template spproj3m*(xx: Spin): untyped =
+  let xp = getPtr xx; template x:untyped {.gensym.} = xp[]
+  let v = [ x[0]-I(x[2]), x[1]+I(x[3]) ]
+  spinVector(2, v)
+template spproj4m*(xx: Spin): untyped =
+  let xp = getPtr xx; template x:untyped {.gensym.} = xp[]
+  let v = [ x[0]-x[2], x[1]-x[3] ]
+  spinVector(2, v)
+
+proc sprecon1p*(r: var auto, x: auto) =
+  ## r: DiracFermion
+  ## x: HalfFermion
+  let nc = x[0].len
+  for i in 0..<nc:
+    r[0][i] = x[0][i]
+    r[1][i] = x[1][i]
+    r[2][i] = x[0][i]
+    r[3][i] = x[1][i]
+
+#[
+template sprecon1p*(x: Spin): untyped = spreconmat1p * x
+template sprecon2p*(x: Spin): untyped = spreconmat2p * x
+template sprecon3p*(x: Spin): untyped = spreconmat3p * x
+template sprecon4p*(x: Spin): untyped = spreconmat4p * x
+template sprecon1m*(x: Spin): untyped = spreconmat1m * x
+template sprecon2m*(x: Spin): untyped = spreconmat2m * x
+template sprecon3m*(x: Spin): untyped = spreconmat3m * x
+template sprecon4m*(x: Spin): untyped = spreconmat4m * x
+]#
+
+template sprecon1p*(xx: Spin): untyped =
+  let xp = getPtr xx; template x:untyped {.gensym.} = xp[]
+  #let v = [ x[0], x[1], -I(x[1]), -I(x[0]) ]
+  #spinVector(4, v)
+  var v {.noInit.}: spinVector(4,evalType(x[0]))
+  #static: echo "sprecon1p: ", v.type
+  v[0] = x[0]
+  v[1] = x[1]
+  v[2] = mI(x[1])
+  v[3] = mI(x[0])
+  v
+template sprecon2p*(xx: Spin): untyped =
+  #let x = xx[]
+  #let v = [ x[0], x[1], x[1], -x[0] ]
+  #spinVector(4, v)
+  let xp = getPtr xx; template x:untyped {.gensym.} = xp[]
+  var v {.noInit.}: spinVector(4,evalType(x[0]))
+  v[0] = x[0]
+  v[1] = x[1]
+  v[2] = x[1]
+  v[3] = -x[0]
+  v
+template sprecon3p*(xx: Spin): untyped =
+  #let x = xx[]
+  #let v = [ x[0], x[1], -I(x[0]), I(x[1]) ]
+  #spinVector(4, v)
+  let xp = getPtr xx; template x:untyped {.gensym.} = xp[]
+  var v {.noInit.}: spinVector(4,evalType(x[0]))
+  v[0] = x[0]
+  v[1] = x[1]
+  v[2] = mI(x[0])
+  v[3] = I(x[1])
+  v
+template sprecon4p*(xx: Spin): untyped =
+  #let x = xx[]
+  #let v = [ x[0], x[1], x[0], x[1] ]
+  #spinVector(4, v)
+  let xp = getPtr xx; template x:untyped {.gensym.} = xp[]
+  var v {.noInit.}: spinVector(4,evalType(x[0]))
+  v[0] = x[0]
+  v[1] = x[1]
+  v[2] = x[0]
+  v[3] = x[1]
+  v
+template sprecon1m*(xx: Spin): untyped =
+  #let x = toRef xx[]
+  #let v = [ x[0], x[1], I(x[1]), I(x[0]) ]
+  #spinVector(4, v)
+  #type T = evalType(x[0])
+  let xp = getPtr xx; template x:untyped {.gensym.} = xp[]
+  var v {.noInit.}: spinVector(4,evalType(x[0]))
+  v[0] = x[0]
+  v[1] = x[1]
+  v[2] = I(x[1])
+  v[3] = I(x[0])
+  v
+template sprecon2m*(xx: Spin): untyped =
+  #let x = xx[]
+  let xp = getPtr xx; template x:untyped {.gensym.} = xp[]
+  let v = [ x[0], x[1], -x[1], x[0] ]
+  spinVector(4, v)
+template sprecon3m*(xx: Spin): untyped =
+  #let x = xx[]
+  let xp = getPtr xx; template x:untyped {.gensym.} = xp[]
+  let v = [ x[0], x[1], I(x[0]), -I(x[1]) ]
+  spinVector(4, v)
+template sprecon4m*(xx: Spin): untyped =
+  #let x = xx[]
+  let xp = getPtr xx; template x:untyped {.gensym.} = xp[]
+  let v = [ x[0], x[1], -x[0], -x[1] ]
+  spinVector(4, v)
+
+when isMainModule:
+  echo gamma0[0,0]
+  let g2 = gamma0 + gamma4
+  #echo g2[0,0]
diff --git a/src/physics/tensorwrap.nim b/src/physics/tensorwrap.nim
index db4fb65..2ed1b50 100644
--- a/src/physics/tensorwrap.nim
+++ b/src/physics/tensorwrap.nim
@@ -17,6 +17,8 @@ template tensorObj*[K,T](k: typedesc[K], t: typedesc[T]): typedesc =
   TensorObj[K,T]
 template tensorObj*[K,T](k: typedesc[K], x: T): auto =
   TensorObj[K,T](obj: x)
+#proc tensorObj*[K,T](k: typedesc[K], x: T): TensorObj[K,getPtr(T)] =
+#  result.obj = getPtr x
 
 template kind*[K,T](x: typedesc[SomeTensor[K,T]]): typedesc = K
 template kind*[K,T](x: SomeTensor[K,T]): typedesc = K
@@ -64,15 +66,30 @@ template `[]`*[K,T,I](x: TensorObj[K,T], i: I): auto =
   else:
     x[][i]
 
-template `[]=`*[K,T,I](x: SomeTensor[K,T], i: I, y: auto) =
+template `[]=`*[K,T,I,Y](x: SomeTensor[K,T], i: I, y: Y) =
+  static: echo "TensorWrap []= KT, I, auto ", $x.type, " ", $I.type, " ", $y.type
   mixin `:=`
-  when I is SomeTensor[K,auto]:
-    x[][i[]] = y
-  elif y is SomeTensor[K,auto]:
+  #when I is SomeTensor2[K,auto]:
+  when I is SomeTensor2:
+    when I.kind is K:
+      x[][i[]] = y
+      return
+  #elif y is SomeTensor2[K,auto]:
+  #elif Y is SomeTensor2 and Y.K is K:
+  when Y is SomeTensor2 and Y.K is K:
     x[][i] = y[]
   else:
     x[][i] = y
 
+#[
+template `[]=`*[K,T,I](x: SomeTensor[K,T], i: I, y: auto) =
+  static: echo "TensorWrap []= KT, I, auto ", $x.type, " ", $I.type, " ", $y.type
+  x[][i] = y
+
+template `[]=`*[K,T,I,Y](x: SomeTensor[K,T], i: I, y: SomeTensor2[K,Y]) =
+  x[][i] = y[]
+]#
+
 template `[]`*[K,T,I](x: TensorObj[K,T], i,j: I): auto =
   x[][i,j]
 template `[]=`*[K,T,I](x: TensorObj[K,T], i,j: I, y: auto) =
@@ -229,16 +246,16 @@ template setBinopTT(f, g) {.dirty.} =
     when X.kind is Y.kind:
       tensorObj(X.kind, f(X[],Y[]))
     elif X.kind.has Y.kind:
-      f(X[], Y)
+      tensorObj(X.kind, f(X[], Y))
     else:  # assume Y.kind has X.kind
-      f(X, Y[])
+      tensorObj(Y.kind, f(X, Y[]))
   template f*[X:SomeTensor,Y:SomeTensor2](x: X, y: Y): auto =
     when X.kind is Y.kind:
       tensorObj(X.kind, f(x[],y[]))
     elif X.kind.has Y.kind:
-      f(x[], y)
+      tensorObj(X.kind, f(x[], y))
     else:  # assume Y.kind has X.kind
-      f(x, y[])
+      tensorObj(Y.kind, f(x, y[]))
   #template f*(x: SomeTensor, y: SomeTensor2): auto =
   #  var r {.noInit.}: f(x.type, y.type)
   #  g(r, x, y)