Exercise solutions of 3dmath

3dmath-exs.ts

const deg2rad = (a: number): number => {
    a%=360;
    a=(a+360)%360;
    return Math.PI*a/180;
};

class mx3 {
    constructor(public rows: [vec3, vec3, vec3]) {}
    transpose(): mx3 {
        return new mx3([
            new vec3(this.rows[0].x,this.rows[1].x,this.rows[2].x),
            new vec3(this.rows[0].y,this.rows[1].y,this.rows[2].y),
            new vec3(this.rows[0].z,this.rows[1].z,this.rows[2].z)
        ]);
    }
    mult(M: mx3): mx3 {
        const MM=M.transpose();
        let rs:vec3[]=[];
        for(let i=0; i<3; ++i) {
            let cs:[number,number,number]=[0,0,0];
            for(let j=0; j<3; ++j) cs[j]=this.rows[i].dotprod(MM.rows[j]);
            rs[i]=vec3.init(cs);
        }
        return new mx3([rs[0],rs[1],rs[2]]);
    }
}

class vec3 {
    constructor(public x: number, public y: number, public z: number) {}
    static init(ns: [number, number, number]) {
        return new vec3(ns[0], ns[1], ns[2]);
    }
    mult(k: number): vec3;
    mult(M: mx3): vec3;
    mult(x: any): vec3 {
        if(typeof x === "number") { return new vec3(this.x*x, this.y*x, this.z*x); }
        return new vec3(
            this.dotprod(x.rows[0]),
            this.dotprod(x.rows[1]),
            this.dotprod(x.rows[2])
        );
    }
    add(that: vec3): vec3 {
        return new vec3(
            this.x+that.x,
            this.y+that.y,
            this.z+that.z
        );
    }
    negate(): vec3 { return this.mult(-1); }
    sub(that: vec3): vec3 { return this.add(that.negate()); }
    dotprod(that: vec3): number { return this.x*that.x+this.y*that.y+this.z*that.z; }
    crossprod(that: vec3): vec3 { return new vec3(this.y*that.z-that.y*this.z, this.x*that.z-that.x*this.z, this.x*that.y-that.x*this.y); }
    abs(): number { return Math.sqrt(this.x*this.x+this.y*this.y+this.z*this.z); }
    parallel(n: vec3): vec3 { return n.mult(this.dotprod(n)); }
    perpendicular(n: vec3): vec3 { return this.sub(this.parallel(n)); }
    // 5.1 rotation
    static rotateM(naxis: vec3, rad: number): mx3 {
        const nx=naxis.x, ny=naxis.y, nz=naxis.z;
        const r0 = new vec3(
            nx*nx*(1-Math.cos(rad))+Math.cos(rad),
            nx*ny*(1-Math.cos(rad))+nz*Math.sin(rad),
            nx*nz*(1-Math.cos(rad))-ny*Math.sin(rad)
        );
        const r1 = new vec3(
            nx*ny*(1-Math.cos(rad))-nz*Math.sin(rad),
            ny*ny*(1-Math.cos(rad))+Math.cos(rad),
            ny*nz*(1-Math.cos(rad))+nx*Math.sin(rad)
        );
        const r2 = new vec3(
            nx*nz*(1-Math.cos(rad))+ny*Math.sin(rad),
            ny*nz*(1-Math.cos(rad))-nx*Math.sin(rad),
            nz*nz*(1-Math.cos(rad))+Math.cos(rad)
        );
        return new mx3([r0,r1,r2]);            
    }
    // 5.2 scale
    static scaleM(naxis: vec3, k: number): mx3 {
        const nx=naxis.x, ny=naxis.y, nz=naxis.z;
        const r0 = new vec3(
            1+(k-1)*nx*nx,
            (k-1)*nx*ny,
            (k-1)*nx*nz
        );
        const r1 = new vec3(
            (k-1)*nx*ny,
            1+(k-1)*ny*ny,
            (k-1)*ny*nz
        );
        const r2 = new vec3(
            (k-1)*nx*nz,
            (k-1)*ny*nz,
            1+(k-1)*nz*nz
        );
        return new mx3([r0,r1,r2]);   
    }
    // 5.3 projection
    static projM(naxis: vec3): mx3 { return vec3.scaleM(naxis, 0); }
    // 5.4 reflection
    static reflM(naxis: vec3): mx3 { return vec3.scaleM(naxis, -1); }
    rotate(naxis: vec3, rad: number): vec3 {
/*      const w = this.crossprod(naxis);
        const vper = this.perpendicular(naxis);
        const vper2 = vper.mult(Math.cos(rad)).add(w.mult(Math.sin(rad)));
        return vper2.add(this.parallel(naxis));*/
        const M = vec3.rotateM(naxis, rad);
        return this.mult(M);
    }
    // 5.2 scale
    scale(naxis: vec3, k: number): vec3 {
        const vper = this.perpendicular(naxis);
        const vpar2 = this.parallel(naxis).mult(k);
        return this.add(naxis.mult(this.dotprod(naxis)).mult(k-1));
    }
}

// Ao->Au->Aw
// Au = Ao_x*p + Ao_y*q + Ao_z*r
// Aw = Au + u_org

const trans_o2u = (Ao: vec3, p: vec3, q: vec3, r: vec3): vec3 => (
    p.mult(Ao.x).add(q.mult(Ao.y)).add(r.mult(Ao.z))
);
const trans_u2w = (uorg: vec3, Au: vec3): vec3 => ( uorg.add(Au) );
const trans_o2w = (Ao: vec3, p: vec3, q: vec3, r: vec3, uorg: vec3): vec3 => (
    trans_u2w(uorg, trans_o2u(Ao,p,q,r))
);

// Aw->Au->Ao
// Au = Aw - u_org
// Ao = [Au∙p, Au∙q, Au∙r] if p,q,r are orthonormal (perpendicular and unit length)

const trans_w2u = (uorg: vec3, Aw: vec3): vec3 => ( Aw.sub(uorg) );
const trans_u2o = (Au: vec3, p: vec3, q: vec3, r: vec3) => (
    vec3.init([Au.dotprod(p), Au.dotprod(q), Au.dotprod(r)])
);
const trans_w2o = (Aw: vec3, p: vec3, q: vec3, r: vec3, uorg: vec3): vec3 => (
    trans_u2o(trans_w2u(uorg, Aw), p, q, r)
);

{
    console.log("Ex 3.6");    
    const Ao: vec3[] = [
        vec3.init([-1,2,0]),
        vec3.init([1,2,0]),
        vec3.init([0,0,0]),
        vec3.init([1,5,0.5]),
        vec3.init([0,5,10])
    ];
    const uorg = vec3.init([1,10,3]);
    const p=vec3.init([0.866,0,-0.5]);
    const q=vec3.init([0,1,0]);
    const r=vec3.init([0.5,0,0.866]);

    const resu = Ao.map((v: vec3) => (trans_o2u(v,p,q,r)));
    const resw = Ao.map((v: vec3) => (trans_o2w(v,p,q,r,uorg)));
    console.log(resu);
    console.log(resw);

    const Aw: vec3[] = [
        vec3.init([1,10,3]),
        vec3.init([0,0,0]),
        vec3.init([2.732,10,2.000]),
        vec3.init([2,11,4]),
        vec3.init([1,20,3])
    ];
    const resuu = Aw.map((v: vec3) => (trans_w2u(uorg,v)));
    const reso = Aw.map((v: vec3) => (trans_w2o(v,p,q,r,uorg)));
    console.log(resuu);
    console.log(reso);
}

{
    console.log("Ex 5.8");
    console.log("Q2",vec3.rotateM(new vec3(1,0,0), deg2rad(-22)));
    console.log("Q3",vec3.rotateM(new vec3(0,1,0), deg2rad(30)));
    console.log("Q4",vec3.rotateM(new vec3(0.267,-0.535,0.802), deg2rad(-15)));
    console.log("Q6",vec3.scaleM(new vec3(0.267,-0.535,0.802), 5));
    console.log("Q7",vec3.projM(new vec3(0.267,-0.535,0.802)));
    console.log("Q8",vec3.reflM(new vec3(0.267,-0.535,0.802)));

    const M1=vec3.rotateM(new vec3(0,1,0), deg2rad(30)),M2=vec3.rotateM(new vec3(1,0,0), deg2rad(-22));
    console.log("Q9-1",M1.mult(M2));
    const M1p=vec3.rotateM(new vec3(0,1,0), deg2rad(-30)),M2p=vec3.rotateM(new vec3(1,0,0), deg2rad(22));
    console.log("Q9-2",M2p.mult(M1p));
}