Tglad
Fractal Molossus
Posts: 703
|
|
« Reply #6 on: December 12, 2009, 12:21:53 AM » |
|
Thanks, here's the code: Vector class: class Vector(){ public: func Init(float xx, float yy, float zz) x = xx y = yy z = zz endfunc float func Dot(Vector other) return x*other.x + y*other.y + z*other.z endfunc func Divide(float f) x = x / f y = y / f z = z / f endfunc func Multiply(Vector a, float f) x = a.x * f y = a.y * f z = a.z * f endfunc func MultiplyEquals(float f) x = x * f y = y * f z = z * f endfunc func Normalise() float mag = sqrt(x*x + y * y + z*z) if mag>0 x = x / mag y = y / mag z = z / mag endif endfunc float func Magnitude() return sqrt(x*x + y*y + z*z) endfunc func Add(Vector a, Vector b) x = a.x + b.x y = a.y + b.y z = a.z + b.z endfunc func AddEquals(Vector a) x = x + a.x y = y + a.y z = z + a.z endfunc func Subtract(Vector a, Vector b) x = a.x - b.x y = a.y - b.y z = a.z - b.z endfunc func Cross(Vector a, Vector b) float _x = (a.y * b.z) - (a.z * b.y) float _y = (a.z * b.x) - (a.x * b.z) float _z = (a.x * b.y) - (a.y * b.x)
x = _x y = _y z = _z endfunc float x float y float z }
Initialisation: Vector corners[4] ; unit length tetrahedron corners corners[0] = new Vector() corners[1] = new Vector() corners[2] = new Vector() corners[3] = new Vector() corners[0].Init(0, 1, 0) corners[1].Init(0, -1/3, sqrt(8/9)) corners[2].Init(-sqrt(2/3), -1/3, -(1/3) * sqrt(2)) corners[3].Init(sqrt(2/3), -1/3, -(1/3) * sqrt(2)) Vector normals[4] Vector crosses[4,3] ; crosses[face index][corner index]. cross product of corners for inside-outside measures Vector C = new Vector() int triangles[4,3] ; triangles[face index, corner of face] triangles[0, 0] = 0 ; wish I knew how to do an initialiser list triangles[0, 1] = 1 triangles[0, 2] = 2 triangles[1, 0] = 0 triangles[1, 1] = 2 triangles[1, 2] = 3 triangles[2, 0] = 0 triangles[2, 1] = 3 triangles[2, 2] = 1 triangles[3, 0] = 2 triangles[3, 1] = 1 triangles[3, 2] = 3 ; for use later Vector point = new Vector() Vector corn0 = new Vector() Vector corn1 = new Vector() Vector corn2 = new Vector() Vector bary = new Vector() int it = 0 repeat normals[it] = new Vector() crosses[it, 0] = new Vector() crosses[it, 1] = new Vector() crosses[it, 2] = new Vector() Vector cr1 = new Vector() cr1.Subtract(corners[triangles[it, 2]], corners[triangles[it, 0]]) Vector cr2 = new Vector() cr2.Subtract(corners[triangles[it, 1]], corners[triangles[it, 0]])
normals[it].Cross(cr1, cr2) normals[it].Normalise() ; leave this divide out below for an alternative where face of tetrahedron has magnitude 1, rather than corners. Produces slightly different shape normals[it].Divide(normals[it].Dot(corners[triangles[it, 0]])) ; normal is not unit length so that a corner point gives dot of 1 crosses[it, 0].Cross(corners[triangles[it, 2]], corners[triangles[it, 1]]) crosses[it, 1].Cross(corners[triangles[it, 0]], corners[triangles[it, 2]]) crosses[it, 2].Cross(corners[triangles[it, 1]], corners[triangles[it, 0]]) until (it=it+1)>=4
int squareMap[4,4,3] ; squareMap[quadrant, face index, corner of quadrant] squareMap[0, 0, 0] = 0 squareMap[0, 0, 1] = 1 squareMap[0, 0, 2] = 2 squareMap[0, 1, 0] = 0 squareMap[0, 1, 1] = 2 squareMap[0, 1, 2] = 3 squareMap[0, 2, 0] = 0 squareMap[0, 2, 1] = 3 squareMap[0, 2, 2] = 1 squareMap[0, 3, 0] = 0 squareMap[0, 3, 1] = 3 squareMap[0, 3, 2] = 1
squareMap[1, 0, 0] = 1 squareMap[1, 0, 1] = 0 squareMap[1, 0, 2] = 3 squareMap[1, 1, 0] = 2 squareMap[1, 1, 1] = 0 squareMap[1, 1, 2] = 1 squareMap[1, 2, 0] = 3 squareMap[1, 2, 1] = 0 squareMap[1, 2, 2] = 2 squareMap[1, 3, 0] = 3 squareMap[1, 3, 1] = 0 squareMap[1, 3, 2] = 2
squareMap[2, 0, 0] = 2 squareMap[2, 0, 1] = 3 squareMap[2, 0, 2] = 0 squareMap[2, 1, 0] = 3 squareMap[2, 1, 1] = 1 squareMap[2, 1, 2] = 0 squareMap[2, 2, 0] = 1 squareMap[2, 2, 1] = 2 squareMap[2, 2, 2] = 0 squareMap[2, 3, 0] = 1 squareMap[2, 3, 1] = 2 squareMap[2, 3, 2] = 0
squareMap[3, 0, 0] = 1 squareMap[3, 0, 1] = 3 squareMap[3, 0, 2] = 2 squareMap[3, 1, 0] = 2 squareMap[3, 1, 1] = 1 squareMap[3, 1, 2] = 3 squareMap[3, 2, 0] = 3 squareMap[3, 2, 1] = 2 squareMap[3, 2, 2] = 1 squareMap[3, 3, 0] = 3 squareMap[3, 3, 1] = 2 squareMap[3, 3, 2] = 1
Inner loop: ; convert c and z to vectors, for convenience C.x = real(cri) C.y = imag(cri) C.z = cj point.x = real(zri) point.y = imag(zri) point.z = zj float magnitude = 1 int ij = 0 ; ij is just a unique name for the iterator over each tetrahedron face repeat magnitude = point.Dot(normals[ij]) ; normals isn't unit length so this gives radial magnitude in shape of tetrahedron if (magnitude > 0) ; on right side of triangle ; TDL super simple way to get barycentric coordinates of triangle bary.x = crosses[ij, 0].Dot(point) bary.y = crosses[ij, 1].Dot(point) bary.z = crosses[ij, 2].Dot(point) if (bary.x >= 0 && bary.y >= 0 && bary.z >= 0) ; correct face is found bary.Divide(bary.x+bary.y+bary.z) ; barycentrics sum to 1 int n = 0 if (bary.x >= 0.5) ; quadrant 0 n = 0 bary.x = bary.x - 0.5 elseif (bary.y >= 0.5) ; quadrant 1 n = 1 bary.y = bary.y - 0.5 elseif (bary.x + bary.y <= 0.5) ; quadrant 2 n = 2 else ; quadrant 3 (middle) n = 3 float baryz = 1 - (bary.x + bary.y) bary.y = 0.5 - bary.x bary.x = 0.5 - baryz endif ; Z^2 + c bary.x = bary.x * 2 ; double the size here bary.y = bary.y * 2 bary.z = 1 - (bary.x + bary.y) magnitude = magnitude*magnitude; corn0.Multiply(corners[squareMap[n, ij, 0]], bary.x) corn1.Multiply(corners[squareMap[n, ij, 1]], bary.y) corn2.Multiply(corners[squareMap[n, ij, 2]], bary.z) point.Add(corn0, corn1) point.AddEquals(corn2) point.MultiplyEquals(magnitude) point.AddEquals(C) ij=10 ; break out of loop since the face was found endif endif until (ij = ij+1)>=4 zri = point.x + flip(point.y) ; convert back to z zj = point.z
|