return function(R)
	require("./3d/geom")(R)
	require("./3d/anim")(R)
	require("./sys/object")(R)
	
	R.poly={}
	R.Tri={}
	R.Model={
		Z_SORT=1,	-- Sort faces by average Z
		CALC_NORMAL=2,	-- Calculate face normals
	}

	-- point index: p1,p2,p3
	-- ex: extra data
	function R.Tri:new(p1,p2,p3,ex,normal)
		local o={
			vx={p1,p2,p3},
			ex=ex,
			normal=normal,
		}
		setmetatable(o,self)
		self.__index=self
		return o
	end

	-- vtxs: V4
	-- tris: tri
	--  p[1]-p[3]
	function R.Model:new()
		local o={
			vtxs={},
			anim={},
			tris={},
			passes={},
		}
		setmetatable(o,self)
		self.__index=self
		return o
	end

	-- v=V4
	-- Returns the index
	function R.Model:addVtx(v)
		local i=#self.vtxs+1
		self.vtxs[i]=v
		return i
	end

	-- tri is a Tri
	-- normal can be unset/nil
	function R.Model:addTri(tri)
		local i=#self.tris+1
		self.tris[i]=tri
		return i
	end

	function R.Model:AddDrawPass(name,tris,priority)
		self.passes[#self.passes+1]={name=name,tris=tris,priority=priority}
		table.sort(self.passes,function(a,b) return a.priority<b.priority end)
	end

	-- Transform the geometry:
	-- (local) >entityMat> (world) >viewMat> (view) >projMat> (screen)
	--	(#TODO: collapse that to viewMat -> projMat coming in to this function?)
	-- Transform the normals by entityMat
	-- anim: nil or {cycle=,time=}
	function R.Model:draw(entityMat,viewMat,projMat,viewport,draw,anim,cbFns)
		cbFns={
			geomPrePassFn=cbFns and cbFns.geomPrePassFn or nil,
			geomPostPassFn=cbFns and cbFns.geomPostPassFn or nil,
		}

		local geomMat=projMat:mulM44C(viewMat)
		geomMat:mulM44(entityMat)

		local vtxs={} -- put the model vtxs in a fresh table so that deformations don't alter source data
		local vtxFrame=nil
		if anim==nil or anim.cycle==nil then
			vtxFrame=self.vtxs
		else
			local animDef=self.anim[anim.cycle]
			local unitTime=math.fmod(anim.time,1)	-- Get the fractional part, 0-1
			local iFrame1,interp=math.modf(unitTime*#animDef.frames)
			local iFrame2=1+(iFrame1+1)%#animDef.frames
			iFrame1=iFrame1+1

			-- Blend frames:
			local frame1vtxs=animDef.frames[iFrame1].vtxs
			local frame2vtxs=animDef.frames[iFrame2].vtxs
			for i=1,#animDef.frames[iFrame1].vtxs do
				vtxs[i]=frame1vtxs[i]:lerp(frame2vtxs[i],interp)
			end
			vtxFrame=vtxs
		end

		if cbFns.geomPrePassFn then
			for i=1,#self.vtxs do
				vtxs[i]=cbFns.geomPrePassFn(i,vtxFrame[i]:clone(),self)
			end
		else
			for i=1,#self.vtxs do
				vtxs[i]=vtxFrame[i]:clone()
			end
		end

		local geom={}
		for i=1,#vtxs do
			local p=vtxs[i]:mulM44C(geomMat,true)

			if cbFns.geomPostPassFn then
				p=cbFns.geomPostPassFn(i,p,self)
			end

			local x,y=viewport:transform(p.e[1],p.e[2])
			geom[i]={x=x,y=y,z=p.e[3]}
		end

		if draw then
			if #self.passes>0 then
				for i=1,#self.passes do
					local pass=self.passes[i]
					local drawFn=draw.passFns[pass.name]
					drawFn=drawFn or R.poly.Draw1Colour(math.random(1,15))
					local draw={
						fn=drawFn,
						flags=0,
					}
					self:drawTris(draw,geom,pass.tris,nil)
				end
			elseif draw.fn then
				self:drawTris(draw,geom,self.tris,entityMat)
			end
		end
	end

	function R.Model:drawTris(draw,geom,tris,entityMat)
		local rTris={}
		for i=1,#tris do
			local t=tris[i]

			-- If this triangle doesn't have a normal then create one...
			-- Should only happen on the first run
			if not t.normal then
				t.normal=R.V4Normal(
					self.vtxs[t.vx[1]],
					self.vtxs[t.vx[2]],
					self.vtxs[t.vx[3]]
				)
			end

			local vScr1=geom[t.vx[1]]
			local vScr2=geom[t.vx[2]]
			local vScr3=geom[t.vx[3]]

			local normal=t.normal
			if draw.flags and draw.flags&R.Model.CALC_NORMAL>0 then
				normal=normal:mulM44RotC(entityMat)
			else
				-- Return a dummy normal pointing towards the viewer
				normal=R.V4:new({0,0,1,0})
			end
			
			rTris[#rTris+1]={
				i=i,
				v={vScr1,vScr2,vScr3},
				normal=normal,
				ex=t.ex,
			}
		end

		if draw.flags and draw.flags&R.Model.Z_SORT>0 then
			-- Should this instead sort by closest Z, rather than average?
			table.sort(rTris,
				function(a,b)
					local az=a.v[1].z+a.v[2].z+a.v[3].z
					local bz=b.v[1].z+b.v[2].z+b.v[3].z
					return az>bz
				end
			)
		end
		
		for i=1,#rTris do
			local t=rTris[i]
			local az=t.v[1].z+t.v[2].z+t.v[3].z
	-- This is broke
	--			if az > 0.0 then
				draw.fn(t)
	--			end
		end
	end

	-- Load
	--	vertices: table (one dimension: x,y,z fixed point values)
	--	triangles: table (one dimension: index of p[1]-p[3])
	function R.Model:load(data)
		local scale=data.scale or 1
		local origin=data.origin or {0,0,0}
		local divide_by=data.divide_by or nil

		if divide_by==nil then
			divide_by=0
			-- Get the largest extent of the model. This will be our bound
			for i=1,#data.vertices do
				local extent=math.max(math.abs(data.vertices[i]))
				divide_by=math.max(divide_by,extent)
			end
		end

		self.vtxs={}
		for i=1,#data.vertices//3 do
			local startindex = (i-1)*3+1
			local x=(data.vertices[startindex])/divide_by-origin[1]
			local y=(data.vertices[startindex+1])/divide_by-origin[2]
			local z=(data.vertices[startindex+2])/divide_by-origin[3]
			local v=R.V4:new({x*scale,y*scale,z*scale,1})
			self:addVtx(v)
		end

		self.anim={}
		if data.vanim then
			for name,frames in pairs(data.vanim) do
				local anim=R.Anim:new(name)
				for i,vertices in ipairs(frames) do
					local frame=R.AnimFrame:new()
					trace(#frames)
					for i=1,#vertices//3 do
						local startindex = (i-1)*3+1
						local x=(vertices[startindex])/divide_by-origin[1]
						local y=(vertices[startindex+1])/divide_by-origin[2]
						local z=(vertices[startindex+2])/divide_by-origin[3]
						local v=R.V4:new({x*scale,y*scale,z*scale,1})
						frame:addVtx(v)
					end
					anim:addFrame(frame)
				end
				self.anim[name]=anim
			end
		end

		self.tris={}
		if data.triangles then
			for i=1,#data.triangles//3 do
				local startindex = (i-1)*3+1
				local p1 = data.triangles[startindex]
				local p2 = data.triangles[startindex+1]
				local p3 = data.triangles[startindex+2]
				local ex={}
				if data.triangles_ex then
					ex=data.triangles_ex[i]
				end
				self:addTri(R.Tri:new(p1,p2,p3,ex,nil))
			end
		end

		if data.drawpasses then
			for name,dp in pairs(data.drawpasses) do
				local tris={}
				for i=1,#dp.triangles//3 do
					local startindex = (i-1)*3+1
					local p1 = dp.triangles[startindex]
					local p2 = dp.triangles[startindex+1]
					local p3 = dp.triangles[startindex+2]
--[[					local ex={}
					if data.triangles_ex then
						ex=data.triangles_ex[i]
					end
--]]
					tris[#tris+1]=R.Tri:new(p1,p2,p3,{},nil)
				end
				self:AddDrawPass(name,tris,dp.priority)
			end
		end
	end

	function R.poly.Draw1Colour(colour)
		return function(tr)
			tri(
				tr.v[1].x,tr.v[1].y,
				tr.v[2].x,tr.v[2].y,
				tr.v[3].x,tr.v[3].y,
				colour
			)
		end
	end

	function R.poly.DrawDebug(tr)
		tri(
			tr.v[1].x,tr.v[1].y,
			tr.v[2].x,tr.v[2].y,
			tr.v[3].x,tr.v[3].y,
			1+(tr.i)%15
		)
	end

	function R.poly.DrawFlat(tr)
		tri(
			tr.v[1].x,tr.v[1].y,
			tr.v[2].x,tr.v[2].y,
			tr.v[3].x,tr.v[3].y,
			tr.ex.c
		)
	end

	function R.poly.DrawWire(colour)
		return function(tr)
			trib(
				tr.v[1].x,tr.v[1].y,
				tr.v[2].x,tr.v[2].y,
				tr.v[3].x,tr.v[3].y,
				colour
			)
		end
	end

	--[[
	-- Render example
	-- Lit normals
	-- (NB normalising light in the 
	--  function each time is inefficient
	--  obvs!)
	local DrawLit=function(t)
		local light=V4:new({0,0,-1,1})
		light:normal()

		local l=1+(7.5+t.normal:dot(light)*7.5)//1
		tri(
			t.v[1].x,t.v[1].y,
			t.v[2].x,t.v[2].y,
			t.v[3].x,t.v[3].y,
			l
		)
	end
	--]]

	function R.poly.DrawLit(tr)

		-- #TODO: Move light and screen out of this function
		-- local light=R.V4:new({
		-- 	math.sin(time()*0.0002*math.pi),0,math.cos(time()*0.0002*math.pi),1
		-- })
		-- light:normal()

		local light=R.V4:new({1,1,0})

		-- assume the camera is looking along the
		-- #TODO: pass the scene's camera vector in here so we can cull properly
		local screen=R.V4:new({0,0,1,0})
		screen:normal()

		-- backface culling(ish), really need to project the normals for that
		-- but it eliminates some of the faces which is handy. Seems ok for the
		-- cube but on our chicken other faces are going missing, making me think
		-- they're possibly not all wound in the same direction?

		-- local dp=tr.normal:dot(screen)
		-- if(dp <= 0) then
		-- 	return
		-- end

		local l=tr.normal:dot(light)*14


		l=math.min(14,math.max(0,l))
		l=l+1
		
		tri(
			tr.v[1].x+8,tr.v[1].y,
			tr.v[2].x+8,tr.v[2].y,
			tr.v[3].x+8,tr.v[3].y,
			l
		)
	end

    function R.poly.DrawTex(u0,v0,u1,v1)
        return function(tr)
            ttri(
                tr.v[1].x,tr.v[1].y,
                tr.v[2].x,tr.v[2].y,
                tr.v[3].x,tr.v[3].y,
                u0+tr.ex.u1*u1,v0+tr.ex.v1*v1,
                u0+tr.ex.u2*u1,v0+tr.ex.v2*v1,
                u0+tr.ex.u3*u1,v0+tr.ex.v3*v1,
                2,
                0,
                -tr.v[1].z,-tr.v[2].z,-tr.v[3].z
            )
        end
    end
end
