> with(plots):with(linalg):
# 
# 
# 
> seq({hlis[j][1],hlis[j][2]},j=1..4);
# Enter a list of polygons and colors which define the object.
> hlis:= [ [[[0,0,0],[0,0,3],[3.5,0,5], [7,0,3],[7,0,0]],`coral`],
> [[[0,0,0],[0,0,3],[0,9,3],[0,9,0]],`brown`],
> [[[7,0,0],[7,0,3],[7,9,3],[7,9,0]],`tan`],
> [[[0,9,0],[0,9,3],[3.5,9,5],[7,9,3],[7,9,0]],`orange`],
> [[[0,0,3],[0,9,3],[3.5,9,5],[3.5,0,5]],`green`],
> [[[7,0,3],[7,9,3],[3.5,9,5],[3.5,0,5]],`cyan`],
> [[[3,9.1,0],[5,9.1,0],[5,9.1,8],[3,9.1,8]],`yellow`],
> [[[3,10.2,0],[5,10.2,0],[5,10.2,8],[3,10.2,8]],`red`],
> [[[3,9.1,0],[3,10.2,0],[3,10.2,8],[3,9.1,8]],`maroon`],
> [[[5,9.1,0],[5,10.2,0],[5,10.2,8],[5,9.1,8]],`magenta`]];
# To check, use MAPLES 3-d  graphics routine.
> n:=nops(hlis);
> display3d([seq( polygonplot3d( hlis[j][1], color=hlis[j][2]),
> j=1..n)],axes=normal,scaling=constrained,projection=0.7);
# The following is a sort routine. It sorts the polygons according to
# their depth in the picture, so that the closest ones to the viewer are
# drawn last. This is the painters algorithm. It does not always produce
# the correct order. We do not take into account intersecting polygons.
# It is also confused when three polygons occur such that the first is
# in front of the second is in front of the third is in front of the
# first.
> sortli:=proc() local  nnn, kk,jj,llt;global kproj; nnn:=nops(kproj);
> for kk from 1 to nnn-1 do; for jj from 2 to nnn do; if
> kproj[jj][3]<kproj[jj-1][3]
> then;llt:=kproj[jj-1];kproj[jj-1]:=kproj[jj];kproj[jj]:=llt;fi;od;od;e
> nd:
# Here is the projection by a general linear transformation given by
# multiplying by the matrix A. G maps the object point "lis" to the
# drawing plane. GP maps the point "lis" to its depth.
> A:=matrix([[2,1,1],[-1,2,1],[1,-1,2]]);det(A);
> G:=lis->[A[1,1]*lis[1]+A[1,2]*lis[2]+A[1,3]*lis[3],
> A[2,1]*lis[1]+A[2,2]*lis[2]+A[2,3]*lis[3],
> A[3,1]*lis[1]+A[3,2]*lis[2]+A[3,3]*lis[3]]:
> GP:=lis->[A[1,1]*lis[1]+A[1,2]*lis[2]+A[1,3]*lis[3],
> A[3,1]*lis[1]+A[3,2]*lis[2]+A[3,3]*lis[3]];
> GH:=lis->A[2,1]*lis[1]+A[2,2]*lis[2]+A[2,3]*lis[3];
> kproj:=[seq([map(GP,hlis[j][1]),hlis[j][2],min(op(map(GH,hlis[j][1])))
> ],j=1..n)]:
> sortli(); display([seq( polygonplot( kproj[j][1], color=kproj[j][2]),
> j=1..n)],axes=normal,scaling=constrained);
# Elevations are plotted here. The front view is gotten by taking the
# width and height, but dropping the depth. 
> GPF:=lis->[lis[1],lis[3]]; GPT:=lis->[lis[1],lis[2]+9];
> GPS:=lis->[10+lis[2],lis[3]]; GPFP:=lis->lis[2];  GPTP:=lis->lis[3]; 
> GPSP:=lis->-lis[1];  nn:=nops(hlis);
> kproj:=[seq([map(GPF,hlis[j][1]),hlis[j][2],min(op(map(GPFP,hlis[j][1]
> )))],j=1..n)]: sortli(); kp1:=kproj:
> kproj:=[seq([map(GPT,hlis[j][1]),hlis[j][2],min(op(map(GPTP,hlis[j][1]
> )))],j=1..n)]: sortli(); kp2:=kproj:
> kproj:=[seq([map(GPS,hlis[j][1]),hlis[j][2],min(op(map(GPSP,hlis[j][1]
> )))],j=1..n)]: sortli(); kp3:=kproj: display([seq( polygonplot(
> kp1[j][1], color=kp1[j][2]), j=1..n), seq( polygonplot( kp2[j][1],
> color=kp2[j][2]), j=1..n), seq( polygonplot( kp3[j][1],
> color=kp3[j][2]), j=1..n),textplot([[3,-1,`Front View`], [15,-1,`Side
> View`], [3,20,`Top View`]])],
> scaling=constrained,title=`Elevations`,axes=none);
# Oblique projection. points that are farther from the viewer are
# plotted correspondingly higher and to the right.
> GPo:=lis->[lis[1]+0.707*lis[2],lis[3]+0.707*lis[2]]; 
> GPoP:=lis->lis[2]-0.001*lis[1];    nn:=nops(hlis);
> kproj:=[seq([map(GPo,hlis[j][1]),hlis[j][2],min(op(map(GPoP,hlis[j][1]
> )))],j=1..n)]: sortli(); display([seq( polygonplot( kproj[j][1],
> color=kproj[j][2]), j=1..n),
> plot([GPo([0,-3,0]),GPo([0,10,0])],color=black)],
> scaling=constrained,title=`Oblique Projection`,axes=normal);
# Special isometric projection. The horizontal variables are drawn 30¡
# above and below the axes.
> GPi:=lis->[0.866*lis[1]+0.866*lis[2],-0.5*lis[1]+0.5*lis[2]+lis[3]]; 
> GPiP:=lis->lis[2]-1.001*lis[1];    nn:=nops(hlis);
> kproj:=[seq([map(GPi,hlis[j][1]),hlis[j][2],min(op(map(GPiP,hlis[j][1]
> )))],j=1..n)]: sortli(); display([seq( polygonplot( kproj[j][1],
> color=kproj[j][2]),
> j=1..n),plot({[GPi([0,0,-3]),GPi([0,0,10])],[GPi([0,-3,0]),GPi([0,10,0
> ])],[GPi([-3,0,0]),GPi([10,0,0])]},color=black)],
> scaling=constrained,title=`Standard Isometric Projection`,axes=none);
> GPi:=lis->[0.866*lis[1]-0.866*lis[2],0.5*lis[1]+0.5*lis[2]+lis[3]]; 
> GPiP:=lis->lis[2]+1.001*lis[1]-lis[3];    nn:=nops(hlis);
> kproj:=[seq([map(GPi,hlis[j][1]),hlis[j][2],min(op(map(GPiP,hlis[j][1]
> )))],j=1..n)]: sortli(); display([seq( polygonplot( kproj[j][1],
> color=kproj[j][2]),
> j=1..n),plot({[GPi([0,0,-3]),GPi([0,0,10])],[GPi([0,-3,0]),GPi([0,10,0
> ])],[GPi([-3,0,0]),GPi([10,0,0])]},color=black)],
> scaling=constrained,title=`Standard Isometric Projection`,axes=none);
# Military projection. the horizontal sections are drawn rotated but
# without distortion. the vertical distances are drawn upwards.
> GPm:=lis->[0.6*lis[1]+0.8*lis[2],-0.8*lis[1]+0.6*lis[2]+lis[3]]; 
> GPmP:=lis->-0.8*lis[1]+0.6*lis[2]-lis[3];    nn:=nops(hlis);
> kproj:=[seq([map(GPm,hlis[j][1]),hlis[j][2],min(op(map(GPmP,hlis[j][1]
> )))],j=1..n)]: sortli(); display([seq( polygonplot( kproj[j][1],
> color=kproj[j][2]),
> j=1..n),plot({[GPm([0,0,-3]),GPm([0,0,10])],[GPm([0,-3,0]),GPm([0,10,0
> ])],[GPm([-3,0,0]),GPm([10,0,0])]},color=black)],
> scaling=constrained,title=`Military Projection`,axes=none);
# The perspective projection is computed.  "ep"  is the eye-point. "cp"
# is the center point in the object toward which the eye is looking. "T"
# is the translation that moves the object so that the eyepoint becomes
# the origin. "R" is the rotation about the vertical axis which moves
# the vector from the new eyepoint=origin to the new center point around
# to the "y-z" axis. "r1" is the horizontal distance eye to center, abd
# "s1" and "c1" are the sine and cosine of this rotation. "S" is the
# rotation about the new "y" axis which moves the origin-new center
# vector to the "y"-axis. "r2" is eye to center distance and "s2" and
# "c2" rare the sine and cosine of this rotation.  The composite
# rotation "SRT=S(R(T( )))" is the rigid motion that takes the object to
# the position such that the eye is at the origin, the eye-center vector
# is along the positive "y" axis, and vertical lines are in the vertical
# "y-z" plane. "GPp" is the perpendicular parallel projection of the
# object to the "x-z" plane. "GPpP" finds the "y" coordinate of the
# point (used for determining which polygons are closest to the eye.
> ep:=[11.0,-15.0,2.0];
> cp:=[3.5,5.0,3.0];lis:='lis';
> T:=lis->[lis[1]-ep[1],lis[2]-ep[2],lis[3]-ep[3]]:
> dp:=T(cp):
> r1:=sqrt(dp[1]^2+dp[2]^2):
> s1:=dp[1]/r1:c1:=dp[2]/r1:
> R:=lis->[c1*lis[1]-s1*lis[2],c1*lis[2]+s1*lis[1],lis[3]]:
> rdp:=R(dp):
> r2:=sqrt(dp[1]^2+dp[2]^2+dp[3]^2):
> s2:=rdp[3]/r2:c2:=rdp[2]/r2:
> S:=lis->[lis[1],-s2*lis[2]+c2*lis[3],c2*lis[2]+s2*lis[3]]:
> S(rdp):
> SRT:=lis->S(R(T(lis))):
> GPp:=lis->[op(1,SRT(lis)),op(2,SRT(lis))]:GPpP:=lis->op(3,SRT(lis)):
> GPp(lis):GPpP(lis):
> nn:=nops(hlis);
> kproj:=[seq([map(GPp,hlis[j][1]),hlis[j][2],min(op(map(GPpP,hlis[j][1]
> )))],j=1..nn)]: sortli(); display([seq( polygonplot( kproj[j][1],
> color=kproj[j][2]),
> j=1..n),plot({[GPp([0,0,-3]),GPp([0,0,10])],[GPp([0,-3,0]),GPp([0,10,0
> ])],[GPp([-3,0,0]),GPp([10,0,0])]},color=black),plot([[0,0]],style=poi
> nt,symbol=circle,color=black)], scaling=constrained,title=`Parallel
> Projection in Eye to Center direction`,axes=none);
# "GPx" is the projective transformation, which computes the "x" and "z"
# coordintes of the point in the object as it would appear on the
# drawing plane "y=1". Farther objects are smaller because "GPx" is
# computed by dividing by the distance.
> GPx:=lis->[op(1,GPp(lis))/GPpP(lis),op(2,GPp(lis))/GPpP(lis)]:
> nn:=nops(hlis):kproj:=[seq([map(GPx,hlis[j][1]),hlis[j][2],min(op(map(
> GPpP,hlis[j][1])))],j=1..nn)]: sortli(); display([seq( polygonplot(
> kproj[j][1], color=kproj[j][2]),
> j=1..n),plot({[GPx([0,0,-3]),GPx([0,0,10])],[GPx([0,-3,0]),GPx([0,10,0
> ])],[GPx([-3,0,0]),GPx([10,0,0])]},color=black),plot([[0,0]],style=poi
> nt,symbol=circle,color=black)], scaling=constrained,title=`Perspective
> Projection in Eye to Center direction`,axes=none);
# Another closer point of view to exaggerate the three point
# perspective.
>  

> ep:=[7.9,-5.0,10.0];
> cp:=[3.5,5.0,4.0];lis:='lis';
> T:=lis->[lis[1]-ep[1],lis[2]-ep[2],lis[3]-ep[3]]:
> dp:=T(cp):
> r1:=sqrt(dp[1]^2+dp[2]^2):
> s1:=dp[1]/r1:c1:=dp[2]/r1:
> R:=lis->[c1*lis[1]-s1*lis[2],c1*lis[2]+s1*lis[1],lis[3]]:
> rdp:=R(dp):
> r2:=sqrt(dp[1]^2+dp[2]^2+dp[3]^2):
> s2:=rdp[3]/r2:c2:=rdp[2]/r2:
> S:=lis->[lis[1],-s2*lis[2]+c2*lis[3],c2*lis[2]+s2*lis[3]]:
> S(rdp):
> SRT:=lis->S(R(T(lis))):
> GPp:=lis->[op(1,SRT(lis)),op(2,SRT(lis))]:GPpP:=lis->op(3,SRT(lis)):
> GPp(lis):GPpP(lis):
> nn:=nops(hlis);
> kproj:=[seq([map(GPp,hlis[j][1]),hlis[j][2],min(op(map(GPpP,hlis[j][1]
> )))],j=1..nn)]: sortli(); display([seq( polygonplot( kproj[j][1],
> color=kproj[j][2]),
> j=1..n),plot({[GPp([0,0,-3]),GPp([0,0,10])],[GPp([0,-3,0]),GPp([0,10,0
> ])],[GPp([-3,0,0]),GPp([10,0,0])]},color=black),plot([[0,0]],style=poi
> nt,symbol=circle,color=black)], scaling=constrained,title=`Parallel
> Projection in Eye to Center direction`,axes=none);
> GPx:=lis->[op(1,GPp(lis))/GPpP(lis),op(2,GPp(lis))/GPpP(lis)]:
> nn:=nops(hlis);
> kproj:=[seq([map(GPx,hlis[j][1]),hlis[j][2],max(op(map(GPpP,hlis[j][1]
> )))],j=1..nn)]: sortli(); display([seq( polygonplot( kproj[j][1],
> color=kproj[j][2]),
> j=1..n),plot({[GPx([0,0,-3]),GPx([0,0,10])],[GPx([0,-3,0]),GPx([0,10,0
> ])],[GPx([-3,0,0]),GPx([10,0,0])]},color=black),plot([[0,0]],style=poi
> nt,symbol=circle,color=black)], scaling=constrained,title=`Perspective
> Projection in Eye to Center direction`,axes=none);
# One point perspective. When the eye point and center point have the
# same "x-z" coordinats, then vertical and horizontal lines in the
# object are parallel to the drawing plane, and the projective
# transformation  maps these lines to parallel horizontal and vertical
# lines.  The perpendicular projection is just the front vies.
> ep:=[6.0,-15.0,2.0];
> cp:=[6.0,5.0,2.0];lis:='lis';
> T:=lis->[lis[1]-ep[1],lis[2]-ep[2],lis[3]-ep[3]]:
> dp:=T(cp):
> r1:=sqrt(dp[1]^2+dp[2]^2):
> s1:=dp[1]/r1:c1:=dp[2]/r1:
> R:=lis->[c1*lis[1]-s1*lis[2],c1*lis[2]+s1*lis[1],lis[3]]:
> rdp:=R(dp):
> r2:=sqrt(dp[1]^2+dp[2]^2+dp[3]^2):
> s2:=rdp[3]/r2:c2:=rdp[2]/r2:
> S:=lis->[lis[1],-s2*lis[2]+c2*lis[3],c2*lis[2]+s2*lis[3]]:
> S(rdp):
> SRT:=lis->S(R(T(lis))):
> GPp:=lis->[op(1,SRT(lis)),op(2,SRT(lis))]:GPpP:=lis->op(3,SRT(lis)):
> GPp(lis):GPpP(lis):
> nn:=nops(hlis);
> kproj:=[seq([map(GPp,hlis[j][1]),hlis[j][2],min(op(map(GPpP,hlis[j][1]
> )))],j=1..nn)]: sortli(); display([seq( polygonplot( kproj[j][1],
> color=kproj[j][2]),
> j=1..n),plot({[GPp([0,0,-3]),GPp([0,0,10])],[GPp([0,-3,0]),GPp([0,10,0
> ])],[GPp([-3,0,0]),GPp([10,0,0])]},color=black),plot([[0,0]],style=poi
> nt,symbol=circle,color=black)], scaling=constrained,title=`Parallel
> Projection in Eye to Center direction`,axes=none);
> GPx:=lis->[op(1,GPp(lis))/GPpP(lis),op(2,GPp(lis))/GPpP(lis)]:
> nn:=nops(hlis);
> kproj:=[seq([map(GPx,hlis[j][1]),hlis[j][2],min(op(map(GPpP,hlis[j][1]
> )))],j=1..nn)]: sortli(); display([seq( polygonplot( kproj[j][1],
> color=kproj[j][2]),
> j=1..n),plot({[GPx([0,0,-3]),GPx([0,0,10])],[GPx([0,-3,0]),GPx([0,10,0
> ])],[GPx([-3,0,0]),GPx([10,0,0])]},color=black),plot([[0,0]],style=poi
> nt,symbol=circle,color=black)], scaling=constrained,title=`One Point
> Perspective `,axes=none);
# Two point perspective. The eye point and the center point are at the
# same height. Thus the vertical lines are rigidly moved to vertical
# lines, and the perspective transformation keeps them vertical and
# parallel. Parallel lines in the "x" direction and in the "y" direction
# are mapped to pencils of lines which converge to their own vanishing
# points. 
> ep:=[16.0,-15.0,2.0];
> cp:=[6.0,5.0,2.0];lis:='lis';
> T:=lis->[lis[1]-ep[1],lis[2]-ep[2],lis[3]-ep[3]]:
> dp:=T(cp):
> r1:=sqrt(dp[1]^2+dp[2]^2):
> s1:=dp[1]/r1:c1:=dp[2]/r1:
> R:=lis->[c1*lis[1]-s1*lis[2],c1*lis[2]+s1*lis[1],lis[3]]:
> rdp:=R(dp):
> r2:=sqrt(dp[1]^2+dp[2]^2+dp[3]^2):
> s2:=rdp[3]/r2:c2:=rdp[2]/r2:
> S:=lis->[lis[1],-s2*lis[2]+c2*lis[3],c2*lis[2]+s2*lis[3]]:
> S(rdp):
> SRT:=lis->S(R(T(lis))):
> GPp:=lis->[op(1,SRT(lis)),op(2,SRT(lis))]:GPpP:=lis->op(3,SRT(lis)):
> GPp(lis):GPpP(lis):
> nn:=nops(hlis);
> kproj:=[seq([map(GPp,hlis[j][1]),hlis[j][2],min(op(map(GPpP,hlis[j][1]
> )))],j=1..nn)]: sortli(); display([seq( polygonplot( kproj[j][1],
> color=kproj[j][2]),
> j=1..n),plot({[GPp([0,0,-3]),GPp([0,0,10])],[GPp([0,-3,0]),GPp([0,10,0
> ])],[GPp([-3,0,0]),GPp([10,0,0])]},color=black),plot([[0,0]],style=poi
> nt,symbol=circle,color=black)], scaling=constrained,title=`Parallel
> Projection in Eye to Center direction`,axes=none);
# 
# 
> GPx:=lis->[op(1,GPp(lis))/GPpP(lis),op(2,GPp(lis))/GPpP(lis)]:
> nn:=nops(hlis);
> kproj:=[seq([map(GPx,hlis[j][1]),hlis[j][2],min(op(map(GPpP,hlis[j][1]
> )))],j=1..nn)]: sortli(); display([seq( polygonplot( kproj[j][1],
> color=kproj[j][2]),
> j=1..n),plot({[GPx([0,0,-3]),GPx([0,0,10])],[GPx([0,-3,0]),GPx([0,10,0
> ])],[GPx([-3,0,0]),GPx([10,0,0])]},color=black),plot([[0,0]],style=poi
> nt,symbol=circle,color=black)], scaling=constrained,title=`Two Point
> Perspective `,axes=none);
> 
>