Spherical Harmonic Shapes

Log in to post a comment.

// You can find the Turtle API reference here: https://turtletoy.net/syntax Canvas.setpenopacity(0.05); // Global code will be evaluated once. const turtle = new Turtle(); const PI = Math.PI; const PI2 = Math.PI*2; const sqrt2 = Math.sqrt(2); const RANDOM = 1; const WIREFRAME = 1 // Change COEF to combine the linear combination of let COEF = [ 0.1, 0.1,0.1,0.1, 1,0.1,0.1,0.1,1, 0.01,0.1,0.1,0.01,0.01,0.01,0.1, 0.1,0.001,0.001,0.001,0.001,0.001,0.001,0.001,1, 0.01,0.01,0.01,0.01,0.001,0.01,0.01,0.01,0.001,0.001,0.001, 1,0.001,0.001,0.001,0.0001,0.001,0.0001,0.001,0.0001,0.001,0.0001,0.001,1, //0.0001,0.0001,0.0001,0.0001,1,0.0001,0.0001,0.0001,0.0001,0.0001,0.0001,0.0001,0.0001,0.0001,0.0001 ]; const LEVEL = Math.sqrt(COEF.length); let scale = [1.0,1.0,1.0] let rotation = quat(0.0,3.1,3.1); let position = [0,0,1.5]; // Begin SH Math function factorial(x) { if(x==0){return 1;} else{return x*factorial(x-1);} } function P(l,m,x) { let pmm = 1.0; if(m>0) { let somx2 = Math.sqrt((1-x)*(1+x)); let fact = 1.0; for(let i=1;i<=m;i++) { pmm *= (-fact)*somx2; fact += 2.0; } } if(l==m) return pmm; let pmmp1 = x*(2*m+1.0)*pmm; if(l==m+1){return pmmp1;} let p11 = 0.0; for(let i = m+2; i<=l;i++) { pll = ((2.0*11-1.0)*x*pmmp1 - (i+m-1)*pmm) / (i-m); pmm = pmmp1; pmmp1 = pll; } return pll; } function K(l,m) { let temp = ((2.0*l+1.0)*factorial(l-m))/(4.0*PI*factorial(l+m)); return Math.sqrt(temp); } function SH(l,m,theta,phi) { // theta in range [0...PI] // phi in range [0...2PI] if(m==0){return K(1,0)*P(l,m,Math.cos(theta));} else if(m>0){return sqrt2*K(l,m)*Math.cos(m*phi)*P(l,m,Math.cos(theta));} else { return sqrt2*K(l,-m)*Math.sin(-m*phi)*P(l,-m,Math.cos(theta)); } } function SHBasisIndex(l,m) { return l*(l+1) + m; } if(RANDOM) { for(let i=0;i<COEF.length;i++) { COEF[i] = COEF[i]*Math.random()*0.5; } } function GetSHScale(theta,phi) { let sum = 0; for(var l=0;l<LEVEL;l++) { for(var m=-l; m<=l;m++) { let index = SHBasisIndex(l,m); sum += COEF[index]*SH(l,m,theta,phi); } } return Math.abs(sum); } // End SH Math // Quaternion from Eular angles // https://en.wikipedia.org/wiki/Conversion_between_quaternions_and_Euler_angles function quat(yaw,pitch,roll) { let cy = Math.cos(yaw*0.5); let sy = Math.sin(yaw*0.5); let cp = Math.cos(pitch*0.5); let sp = Math.sin(pitch*0.5); let cr = Math.cos(roll*0.5); let sr = Math.sin(roll*0.5); let w = cy * cp * cr + sy * sp * sr; let x = cy * cp * sr - sy * sp * cr; let y = sy * cp * sr + cy * sp * cr; let z = sy * cp * cr - cy * sp * sr; return [x,y,z,w]; } // ref https://github.com/toji/gl-matrix/blob/master/src/mat4.js // fromRotationTranslationScale // matrix is stored in the row pattern function make_transform(o,pos,quat,scale) { // scale temp vars let sx = scale[0]; let sy = scale[1]; let sz = scale[2]; // translation temp vars let tx = pos[0]; let ty = pos[1]; let tz = pos[2]; //Quat calc temp variables let x = quat[0]; let y = quat[1]; let z = quat[2]; let w = quat[3]; let xx = x*x*2.0; let xy = x*y*2.0; let xz = x*z*2.0; let wx = x*w*2.0; let yy = y*y*2.0; let yz = y*z*2.0; let wy = y*w*2.0; let zz = z*z*2.0; let wz = z*w*2.0; o[0] = (1-yy-zz)*sx; o[1] = (xy-wz)*sy; o[2] = (xz+wy)*sz; o[3] = tx; o[4] = (xy+wz)*sx; o[5] = (1-xx-zz)*sy; o[6] = (yz-wx)*sz; o[7] = ty; o[8] = (xz-wy)*sx; o[9] = (yz+wx)*sy; o[10] = (1-xx-yy)*sz; o[11] = tz; o[12] = 0; o[13] = 0; o[14] = 0; o[15] = 1; return o; } function transform_vector(v,m) { let x = v[0], y = v[1], z = v[2], w = v[3]; let o = [0, 0, 0, 0]; o[0] = m[0] * x + m[1] * y + m[2] * z + m[3] * w; o[1] = m[4] * x + m[5] * y + m[6] * z + m[7] * w; o[2] = m[8] * x + m[9] * y + m[10] * z + m[11] * w; o[3] = m[12] * x + m[13] * y + m[14] * z + m[15] * w; return o; } let focus = 200; function project(v){ let x = v[0]; let y = v[1]; let z = v[2]; let s = focus/z; return [x*s,y*s]; } //transfrom matrix let transform = new Float32Array(16); let Mesh = function(verts,edges,tsf) { let out = {} out.verts = verts; out.edges = edges; out.transform = tsf return out; } // fill in the SH shape vertice and edge list function make_sh_shape(nu,nv) { let o = {}; o.edges = []; o.verts = []; let du = 1.0/nu; let dv = 1.0/nv; nu +=1; nv +=1; let u = 0,v = 0; let i = 0,j = 0; if(!WIREFRAME) { o.verts.push([0,0,0,1.0]); } // push vertice for(let i=0;i<nu;i++) { u = i*du; for(let j=0;j<nv;j++) { v = j*dv; let cu = Math.cos(u*PI2); let su = Math.sin(u*PI2); let cv = Math.cos(v*PI); let sv = Math.sin(v*PI); let x = sv*cu; let y = sv*su; let z = cv; let sh_scale = GetSHScale(v*PI,u*PI2); o.verts.push([x*sh_scale,y*sh_scale,z*sh_scale,1.0]); //o.verts.push([u,v,1.0,1.0]); } } //link edges for(let i=0;i<nu;i++) { u = i*du; for(let j=0;j<nv;j++) { v = j*dv; let id = i*nv+j; let inext = (i+1)%nu; let jnext = (j+1)%nv; let id_down = i*nv+jnext; let id_right = inext*nv+j; if(WIREFRAME) { o.edges.push(id,id_down); o.edges.push(id,id_right); } else { o.edges.push(0,id); } } } return o; } function draw_mesh(tctx,mesh) { let edges = mesh.edges; let verts = mesh.verts; let tsf = mesh.transform; for (let i = 0, len = edges.length/2; i < len; i++) { let i0 = edges[2*i]; let i1 = edges[2*i+1]; // perform transform let v0 = transform_vector(verts[i0],tsf); let v1 = transform_vector(verts[i1],tsf); // projection v0 = project(v0); v1 = project(v1); // draw edge tctx.penup(); tctx.goto([v0[0],v0[1]]); tctx.pendown(); tctx.goto([v1[0],v1[1]]); } } function draw_mesh_substep(tctx,mesh,step) { let edges = mesh.edges; let len = edges.length/2; if(step>=len-2) { return false; } let verts = mesh.verts; let tsf = mesh.transform; let i0 = edges[2*step]; let i1 = edges[2*step+1]; // perform transform let v0 = transform_vector(verts[i0],tsf); let v1 = transform_vector(verts[i1],tsf); // projection v0 = project(v0); v1 = project(v1); // draw edge tctx.penup(); tctx.goto([v0[0],v0[1]]); tctx.pendown(); tctx.goto([v1[0],v1[1]]); return true; } make_transform(transform,position,rotation,scale); let btl = make_sh_shape(200,100); let bottle = Mesh(btl.verts,btl.edges,transform); // The walk function will be called until it returns false. function walk(i) { //draw_mesh(turtle,bottle); let stop = draw_mesh_substep(turtle,bottle,i); return stop; }