summaryrefslogtreecommitdiff
path: root/Master/texmf-dist/metapost/featpost
diff options
context:
space:
mode:
authorKarl Berry <karl@freefriends.org>2014-10-10 22:04:00 +0000
committerKarl Berry <karl@freefriends.org>2014-10-10 22:04:00 +0000
commit2afe81e876627e6ab68d6f8999467e2559ad70b1 (patch)
tree94a7ba24efbcc65ac990ad2e2a7a4c5790fce4a8 /Master/texmf-dist/metapost/featpost
parentca3d251e96e4ab3ac745a18d1a8ef843a8cd803a (diff)
featpost (10oct14)
git-svn-id: svn://tug.org/texlive/trunk@35346 c570f23f-e606-0410-a88d-b1316a301751
Diffstat (limited to 'Master/texmf-dist/metapost/featpost')
-rw-r--r--Master/texmf-dist/metapost/featpost/featpost3Dplus2D.mp327
1 files changed, 249 insertions, 78 deletions
diff --git a/Master/texmf-dist/metapost/featpost/featpost3Dplus2D.mp b/Master/texmf-dist/metapost/featpost/featpost3Dplus2D.mp
index 8a005bfd6d4..52a5f7789ed 100644
--- a/Master/texmf-dist/metapost/featpost/featpost3Dplus2D.mp
+++ b/Master/texmf-dist/metapost/featpost/featpost3Dplus2D.mp
@@ -7,7 +7,7 @@
% P. Jørgensen
% S. Pakin
%
-% Copyright (C) 2013
+% Copyright (C) 2014
% This set of macros adds a lot of features to
% the MetaPost language and eases the production of
@@ -23,7 +23,7 @@
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
% GNU General Public License for more details.
- message "Preloading FeatPost macros, version 0.8.7";
+ message "FeatPost-0.8.8";
warningcheck := 0;
background := 0.987white;
@@ -204,7 +204,6 @@
enddef;
% The normalized crossproduct of two vectors.
-% Also check getangle below.
def ncrossprod(expr A, B) =
N( ccrossprod( A, B ) )
@@ -290,6 +289,18 @@
endgroup
enddef;
+% Maybe you would like to calculate the angular arguments of kindofcube...
+
+ def getanglepair( expr InVec ) =
+ begingroup
+ save alphaone, alphatwo;
+ numeric alphaone, alphatwo;
+ alphaone = angle( ( X(InVec), Y(InVec) ) );
+ alphatwo = angle( ( X(InVec) ++ Y(InVec), Z(InVec) ) );
+ ( (alphaone,alphatwo) )
+ endgroup
+ enddef;
+
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%% Auxiliary:
@@ -366,6 +377,8 @@
endgroup
enddef;
+% The following two procedures are useful for getready.
+
vardef cstr( expr Cl ) =
"(" &
decimal(X(Cl)) &
@@ -436,7 +449,8 @@
enddef;
% Much improved rigorous pseudo-projection algorithm that follows
-% an idea from Cristian Barbarosie. This makes shadows.
+% an idea from Cristian Barbarosie.
+% This makes shadows caused by a light source point.
def cb(expr R) =
begingroup
@@ -477,6 +491,14 @@
endgroup
enddef;
+% Vanishing point.
+
+ def vp( expr D ) =
+ begingroup
+ ( rp( lineintersectplan( f, D, viewcentr, f-viewcentr) ) )
+ endgroup
+ enddef;
+
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%% Basic Functions:
@@ -1426,6 +1448,10 @@
% And again wrong. The last is tdcircarrow.
% Well, what can I say, really the last is ellipsoid.
% Wait! It is torushadow!
+% Bullshit. Only death can stop me. Now it's revolparab.
+% Continuing with torus accessories...
+% And I forgot to mention intersectprolatespheroid!!!
+% The first of Red October, 2014, added vp.
def spheroidshadow( expr CentrPoi, NorthPoleVec, Ray ) =
begingroup
@@ -1528,7 +1554,8 @@
endgroup
enddef;
-% Another brute-force algorythm. It's advisable to use three orthogonal axes.
+% Another brute-force algorythm.
+% It's advisable to use three orthogonal axes.
def ellipsoid( expr Centr, AxOne, AxTwo, AxThr ) =
begingroup
@@ -1600,6 +1627,104 @@
endgroup
enddef;
+ def revolparab( expr BaseCenter, ParabTip, BaseRay ) =
+ begingroup
+ save bcpt, conetip, fakex, fakey, fakez, tipview, ellicenter, coneview;
+ save tanefe, majorvec, minorvec, cutvec, auxpoi, crux, auxx, baseview;
+ save xzero, yzero, xdelta, fakea, xpos, ypos, xneg, yneg, l, apertur;
+ save auxy, maxy, ellmaxang, ymin, xefe, auxray, auxcos, auxsin, ste, a;
+ save auxpath, tippath;
+ save conda, condb, condc;
+ color bcpt, conetip, fakex, fakey, fakez, tipview, ellicenter, coneview;
+ color tanefe, majorvec, minorvec, cutvec, auxpoi, crux, auxx, baseview;
+ numeric xzero, yzero, xdelta, fakea, xpos, ypos, xneg, yneg, l, apertur;
+ numeric auxy, maxy, ellmaxang, ymin, xefe, auxray, auxcos, auxsin, ste;
+ numeric a;
+ path auxpath, tippath;
+ boolean conda, condb, condc;
+ bcpt = BaseCenter-ParabTip;
+ conetip = BaseCenter-2*bcpt;
+ maxy = conorm(bcpt);
+ apertur = angle(2*maxy,BaseRay);
+ fakey = N(bcpt);
+ coneview = f-conetip;
+ tipview = f-ParabTip;
+ baseview = f-BaseCenter;
+ a = getangle(coneview,bcpt);
+ conda = a<=apertur;
+ condb = a>=180-apertur;
+ fakez = ncrossprod( tipview, bcpt );
+ fakex = ccrossprod( fakey, fakez );
+ xzero = cdotprod( fakex, tipview );
+ yzero = cdotprod( fakey, tipview );
+ fakea = maxy/(BaseRay**2);
+ condc = yzero>=fakea*(xzero**2);
+ if (conda and condc) or condb:
+ auxpath = rigorouscircle( BaseCenter, bcpt, BaseRay );
+ unfill auxpath;
+ draw auxpath;
+ else:
+ xdelta = sqrt(xzero**2 - yzero/fakea);
+ xpos = xzero + xdelta;
+ xneg = xzero - xdelta;
+ ypos = fakea*(xpos**2);
+ yneg = fakea*(xneg**2);
+ auxy = 0.5[ypos,yneg];
+ ellicenter = ParabTip+fakex*xzero+fakey*auxy;
+ if yneg<ypos:
+ ymin=yneg;
+ xefe=xneg;
+ else:
+ ymin=ypos;
+ xefe=xpos;
+ fi;
+ tanefe = ParabTip +fakex*xefe+fakey*ymin;
+ minorvec = xdelta*fakez;
+ majorvec = tanefe-ellicenter;
+ if (yneg<maxy) and (ypos<maxy):
+ draw ellipticpath( ellicenter, majorvec, minorvec );
+ else:
+ auxpoi = lineintersectplan(conetip,conetip-f,BaseCenter,bcpt);
+ auxray = conorm( BaseCenter-auxpoi );
+ auxcos = BaseRay/auxray;
+ auxsin = 1 +-+ auxcos;
+ if cdotprod(fakex,auxpoi-BaseCenter)>0:
+ auxx = fakex;
+ else:
+ auxx = -fakex;
+ fi;
+ crux = BaseCenter+BaseRay*(auxx*auxcos +fakez*auxsin);
+ cutvec = crux-ellicenter;
+ auxcos := cdotprod(cutvec,N(majorvec))/conorm(majorvec);
+ auxsin := cdotprod(cutvec,N(minorvec))/conorm(minorvec);
+ ellmaxang = angle(auxcos,auxsin);
+ ste = ellmaxang/18; %%%%%%%%%%%%%%%%%%%%%%%%%%%% DANGER %%%%%%%%%%%%
+ tippath =
+ rp(ellicenter+planarrotation(majorvec,minorvec,-ellmaxang))
+ for l=ste-ellmaxang step ste until ellmaxang:
+ ..rp(ellicenter+planarrotation(majorvec,minorvec,l))
+ endfor --rp(ellicenter+planarrotation(majorvec,minorvec,ellmaxang));
+ if cdotprod( baseview, bcpt )<0:
+ auxpath = rigorouscone(true,BaseCenter,bcpt,BaseRay,conetip);
+ numeric len;
+ len = length auxpath;
+ auxpath := subpath (1,len-1) of auxpath;
+ auxpath := tippath--(reverse auxpath)--cycle;
+ unfill auxpath;
+ draw auxpath;
+ else:
+ tippath := tippath--cycle;
+ unfill tippath;
+ draw tippath;
+ auxpath = rigorouscircle( BaseCenter, bcpt, BaseRay );
+ unfill auxpath;
+ draw auxpath;
+ fi;
+ fi;
+ fi
+ endgroup
+ enddef;
+
% You can't see through this hole. f must not be on the hole axis.
% Not yet documented because "buildcycle" doesn't work properly.
@@ -1748,18 +1873,6 @@
endgroup
enddef;
-% Maybe you would like to calculate the angular arguments of kindofcube...
-
- def getanglepair( expr InVec ) =
- begingroup
- save alphaone, alphatwo;
- numeric alphaone, alphatwo;
- alphaone = angle( ( X(InVec), Y(InVec) ) );
- alphatwo = angle( ( X(InVec) ++ Y(InVec), Z(InVec) ) );
- ( (alphaone,alphatwo) )
- endgroup
- enddef;
-
% It's a bit late now but the stage must be set.
def setthestage( expr NumberOfSideSquares, SideSize ) =
@@ -1855,6 +1968,9 @@
endgroup
enddef;
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%%%% Toroidal Stuff
+
% Take a donut.
def smoothtorus( expr Tcenter, Tmoment, Bray, Sray ) =
@@ -2077,6 +2193,54 @@
endgroup
enddef;
+ def pointinsidetorus( expr Point, Tcenter, Tmoment, Bray, Sray ) =
+ begingroup
+ save outputboolean, height, dist, dray;
+ boolean outputboolean;
+ numeric height, dist, dray;
+ height = cdotprod(N(Tmoment),Point-Tcenter);
+ if abs(height)>=Sray:
+ outputboolean = false;
+ else:
+ dist = conorm(Point-Tcenter-height*N(Tmoment));
+ dray = Sray +-+ abs(height);
+ if (dist<Bray+dray) and (dist>Bray-dray):
+ outputboolean = true;
+ else:
+ outputboolean = false;
+ fi;
+ fi;
+ ( outputboolean )
+ endgroup
+ enddef;
+
+ def pointrelativetotorus( expr Point, Tcenter, Tmoment, Bray, Sray ) =
+ begingroup
+ save height, dist;
+ numeric height;
+ color dist;
+ height = cdotprod(N(Tmoment),Point-Tcenter);
+ dist = N(Point-Tcenter-height*N(Tmoment));
+ ( conorm(Point-Bray*dist)-Sray )
+ endgroup
+ enddef;
+
+ def intersectorus( expr Tcenter, Tmoment, Bray, Sray, LinePoi, LineDir ) =
+ begingroup
+ save trypoi, factry, linedi;
+ color trypoi, linedi;
+ numeric factry, j, auxd;
+ trypoi = LinePoi;
+ factry = 0.25;
+ linedi = N(LineDir);
+ for j= 1 upto 50:
+ auxd := pointrelativetotorus( trypoi, Tcenter, Tmoment, Bray, Sray );
+ trypoi := trypoi+factry*linedi*auxd;
+ endfor;
+ ( trypoi )
+ endgroup
+ enddef;
+
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%% Non-standard objects:
@@ -2257,7 +2421,7 @@
enddef;
% Another point is that I want to draw trajectories in space and
-% dependant on velocity: VecFunc( position, velocity ).
+% dependent on velocity: VecFunc( position, velocity ).
% This time is fourth-order Runge-Kutta.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% CHANGES PrintStep!!!!
@@ -3537,7 +3701,7 @@
enddef;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-% Part IV (automatic perspective tuning, polyhedric vertex approximation):
+% Part IV (automatic perspective tuning and minimization):
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
def randomfear =
@@ -3690,6 +3854,67 @@
endgroup
enddef;
+% Minimization routine for scalar functions like y=f(x) where an initial
+% triplet (x1,x2,x3) with x1<x2<x3 is given as a parabolic squeleton that
+% provides a way to search for the smallest value of y (if iterated)
+
+ def minimizestep( expr Abcisscolor )( text PlainFunc ) =
+ begingroup
+ save xa, xb, xc, xd, ya, yb, yc, yd, aux, coeb, coec, den;
+ save colout;
+ numeric xa, xb, xc, xd, ya, yb, yc, yd, aux, coeb, coec, den;
+ color colout;
+ xa = X( Abcisscolor );
+ xb = Y( Abcisscolor );
+ xc = Z( Abcisscolor );
+ ya = PlainFunc(xa);
+ yb = PlainFunc(xb);
+ yc = PlainFunc(xc);
+ if ya = yb:
+ colout = (-0.125[xa,xb],xb,xc);
+ elseif yb = yc:
+ colout = (xa,xb,1.125[xb,xc]);
+ else:
+ if (yb>ya) or (yb>yc):
+ show Abcisscolor;
+ message " Unable to minimizestep!";
+ fi;
+ den = (xb-xc)*((xa**2)-(xb**2))-(xa-xb)*((xb**2)-(xc**2));
+ if abs(den) < 0.0005:
+ show den;
+ message " Unable to minimizestep!";
+ fi;
+ coeb = ((yb-yc)*((xa**2)-(xb**2))-(ya-yb)*((xb**2)-(xc**2)))/den;
+ coec = ((xb-xc)*(ya-yb)-(xa-xb)*(yb-yc))/den;
+ xd = -0.5*coeb/coec;
+ yd = PlainFunc( xd );
+ if ((xa<xd) and (xd<xb)):
+ if (yd<yb):
+ colout = (xa,xd,xb);
+ else:
+ colout = (xd,xb,xc);
+ fi;
+ elseif ((xb<xd) and (xd<xc)):
+ if (yd<yb):
+ colout = (xb,xd,xc);
+ else:
+ colout = (xa,xb,xd);
+ fi;
+ else:
+ aux := 0.125[xb,xc]-0.125[xb,xa];
+ colout = (xa,0.125[xb,xa]+uniformdeviate(aux),xc);
+ fi;
+ fi;
+ ( colout )
+ endgroup
+ enddef;
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%%%% Part V:
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%%%% Intersections:
+%%%%%%%%%%%%%%%%%%%%%%%%%%%% Check lineintersectplan above %%%%%%%%%%%%%%%
+
def calculatecostver(expr VerA,DisA,VerB,DisB,VerC,DisC,TryV) =
begingroup
save a, b, c;
@@ -3700,7 +3925,8 @@
( ( a ++ b ++ c )**2 )
endgroup
enddef;
-
+
+% Approximation of the intersection of three spheres.
% Be aware of the next three danger parameters.
def improvertex( expr VerA, DisA, VerB, DisB, VerC, DisC, IniV ) =
@@ -3861,64 +4087,10 @@
endgroup
enddef;
-% Minimization routine for scalar functions like y=f(x) where an initial
-% triplet (x1,x2,x3) with x1<x2<x3 is given as a parabolic squeleton that
-% provides a way to search for the smallest value of y (if iterated)
-
- def minimizestep( expr Abcisscolor )( text PlainFunc ) =
- begingroup
- save xa, xb, xc, xd, ya, yb, yc, yd, aux, coeb, coec, den;
- save colout;
- numeric xa, xb, xc, xd, ya, yb, yc, yd, aux, coeb, coec, den;
- color colout;
- xa = X( Abcisscolor );
- xb = Y( Abcisscolor );
- xc = Z( Abcisscolor );
- ya = PlainFunc(xa);
- yb = PlainFunc(xb);
- yc = PlainFunc(xc);
- if ya = yb:
- colout = (-0.125[xa,xb],xb,xc);
- elseif yb = yc:
- colout = (xa,xb,1.125[xb,xc]);
- else:
- if (yb>ya) or (yb>yc):
- show Abcisscolor;
- message " Unable to minimizestep!";
- fi;
- den = (xb-xc)*((xa**2)-(xb**2))-(xa-xb)*((xb**2)-(xc**2));
- if abs(den) < 0.0005:
- show den;
- message " Unable to minimizestep!";
- fi;
- coeb = ((yb-yc)*((xa**2)-(xb**2))-(ya-yb)*((xb**2)-(xc**2)))/den;
- coec = ((xb-xc)*(ya-yb)-(xa-xb)*(yb-yc))/den;
- xd = -0.5*coeb/coec;
- yd = PlainFunc( xd );
- if ((xa<xd) and (xd<xb)):
- if (yd<yb):
- colout = (xa,xd,xb);
- else:
- colout = (xd,xb,xc);
- fi;
- elseif ((xb<xd) and (xd<xc)):
- if (yd<yb):
- colout = (xb,xd,xc);
- else:
- colout = (xa,xb,xd);
- fi;
- else:
- aux := 0.125[xb,xc]-0.125[xb,xa];
- colout = (xa,0.125[xb,xa]+uniformdeviate(aux),xc);
- fi;
- fi;
- ( colout )
- endgroup
- enddef;
-
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-%%%% Part V (strictly two-dimensional):
+%%%% Part VI (strictly two-dimensional):
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
% Verify if a path is cyclic (written by Scott Pakin)
def is_cyclic expr cpath =
@@ -4246,8 +4418,7 @@
enddef;
% Shrink or swell a cyclic path without cusp points and without
-% coinciding pre and post control points. This algorithm should
-% be improved to add circular arcs on the outside of convex corners.
+% coinciding pre and post control points.
def lasermachine( expr DefinedPath, Beam, CosLimit ) =
begingroup