Download Learning mEtapost by Doing

\begin{document}
etex
% u: dimensional unit
% xp, yp, zp: coordinates of light source
% bf : brightness factor
% base_color : base color
numeric u,xp,yp,zp,bf;
color base_color;
u = 1cm;
xp := 3; yp := 3; zp := 5;
bf := 30;
base_color := red+green;
% O, Xr, Yr, Zr : reference frame
pair O,Xr,Yr,Zr;
O = (0,0);
Xr = (-.7,-.7) scaled u;
Yr = (1,0) scaled u;
Zr = (0,1) scaled u;
% for drawing the reference frame
vardef frameXYZ(expr s) =
drawarrow O--Xr scaled s;
drawarrow O--Yr scaled s;
drawarrow O--Zr scaled s;
label.llft(btex $X$ etex scaled 1.25, (Xr scaled s));
label.rt(btex $Y$ etex scaled 1.25, (Yr scaled s));
label.top(btex $Z$ etex scaled 1.25, (Zr scaled s));
enddef;
% from 3D to 2D coordinates
vardef project(expr x,y,z) = x*Xr + y*Yr + z*Zr enddef;
% numerical derivatives by central differences
vardef diffx(suffix f)(expr x,y) =
numeric h; h := 0.01;
(f(x+h,y)-f(x-h,y))/(2*h)
enddef;
vardef diffy(suffix f)(expr x,y) =
numeric h; h := 0.01;
(f(x,y+h)-f(x,y-h))/(2*h)
enddef;
% Compute brightness factor at a point
vardef brightnessfactor(suffix f)(expr x,y,z) =
numeric dfx,dfy,ca,cb,cc;
69