summaryrefslogtreecommitdiff
path: root/Master/texmf-dist/metapost
diff options
context:
space:
mode:
Diffstat (limited to 'Master/texmf-dist/metapost')
-rw-r--r--Master/texmf-dist/metapost/featpost/featpost.mp36
-rw-r--r--Master/texmf-dist/metapost/featpost/featpost3Dplus2D.mp1525
2 files changed, 1132 insertions, 429 deletions
diff --git a/Master/texmf-dist/metapost/featpost/featpost.mp b/Master/texmf-dist/metapost/featpost/featpost.mp
deleted file mode 100644
index 9c30b765fed..00000000000
--- a/Master/texmf-dist/metapost/featpost/featpost.mp
+++ /dev/null
@@ -1,36 +0,0 @@
-% featpost.mp
-% L. Nobre G.
-% 2005
-%
-% Use "export TEX=latex" in your .bashrc.
-% Begin by pre-compiling this set of macros with "inimpost featpost.mp".
-% Produce your PostScript (PS) figures with "mpost -mem featpost file".
-% The output of this command (one or several files named "file.N")
-% must be in the current directory to use the bashscripts laproof,
-% lbproof and lcproof. Produce your EncapsulatedPostScript (EPS) figures
-% with "bashscript/lXproof file N". laproof tranforms PS into EPS if
-% the figure is smaller than an A4 page. lbproof produces JPEG and EPS
-% that fits the width of an A4 portrait page and lcproof produces EPS
-% that fits the height of an A4 portrait page independently of its
-% original size.
-% N is the number of the figure (in file) that you want to encapsulate.
-% It is not necessary to use "export TEX=latex" nor "lXproof file N"
-% if the figures have no text.
-% The resulting EPS figures are not insertable in LaTex documents.
-% The originals file.N are.
-
-input plain;
-input featpost3Dplus2D;
-
-dump;
-
-% It is possible to interactively experiment the effect
-% of figure parameters. This requires:
-% 1) MetaPost code written as "anglinerigorouscircle.mp".
-% 2) xcmd package (kindly provided by Pedro Sebastião);
-% it may be downloaded from http://lince.cii.fc.ul.pt/
-% in Debian, Red Hat or tar.gz package formats.
-% 3) Perl.
-% 4) gv.
-% 5) run command "xcmd/xmpost anglinerigorouscircle" or
-% "xcmd/xmpost anymetapostwoextension".
diff --git a/Master/texmf-dist/metapost/featpost/featpost3Dplus2D.mp b/Master/texmf-dist/metapost/featpost/featpost3Dplus2D.mp
index b6055afae58..3c1c809eabd 100644
--- a/Master/texmf-dist/metapost/featpost/featpost3Dplus2D.mp
+++ b/Master/texmf-dist/metapost/featpost/featpost3Dplus2D.mp
@@ -1,17 +1,21 @@
% featpost3Dplus2D.mp
-% L. Nobre G., C. Barbarosie, J. Schwaiger and B. Jackowski
-% nobre@lince.cii.fc.ul.pt
-% http://matagalatlante.org
-% Copyright (C) 2005
-% see also featpost.mp
-
-% This set of macros extends the MetaPost language
-% to three dimensions and eases the production of
-% physics diagrams.
+% L. Nobre G., lnobreg@gmail.com, http://matagalatlante.org
+% C. Barbarosie
+% J. Schwaiger
+% B. Jackowski
+% P. J. Sebastião
+% P. Jørgensen
+% S. Pakin
+%
+% Copyright (C) 2011
+
+% This set of macros adds a lot of features to
+% the MetaPost language and eases the production of
+% physics diagrams and animations.
% This is free software; you can redistribute it and/or
% modify it under the terms of the GNU General Public License
-% as published by the Free Software Foundation; either version 2
+% as published by the Free Software Foundation; either version 3
% of the License, or (at your option) any later version.
% This set of macros is distributed in the hope that it will be useful,
@@ -19,16 +23,17 @@
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
% GNU General Public License for more details.
- message "Preloading FeatPost macros, version 0.6.7";
+ message "Preloading FeatPost macros, version 0.8.2";
warningcheck := 0;
-
background := 0.987white;
+ defaultscale := 0.75;
+ defaultfont := "cmss17"; % This is used by cartaxes
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% Global Variables %%%%%%%%%%%
boolean ParallelProj, SphericalDistortion, FCD[], ShadowOn;
- boolean OverRidePolyhedricColor;
+ boolean OverRidePolyhedricColor, MalcomX;
numeric Nobjects, RefDist[], HoriZon, RopeColorSeq[], PhotoMarks;
numeric Spread, PrintStep, PageHeight, PageWidth, ActuC, Shifts;
numeric NL, npl[], NF, npf[], FC[], MaxFearLimit, TableColors;
@@ -40,10 +45,6 @@
string ostr[];
pen BackPen, ForePen;
-%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% Kept for backward compatibility
-
- Shifts := 105.00mm;
-
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% Default Values %%%%%%%%%%%%%%%
f := (3,5,4); % This f is the point of view in 3D
@@ -52,10 +53,12 @@
Spread := 140; % Magnification
- ShiftV := 105.00mm*(1,1); % Central coordinates on paper
+ Shifts := 105.00mm;
+
+ ShiftV := Shifts*(1,1); % Central coordinates on paper
- OriginProjPagePos := (105.00mm,148.45mm); % This should be the
- % page center.
+ OriginProjPagePos := (Shifts,148.45mm); % This should be the
+ % page center.
ParallelProj := false; % Kind of perspective
% Can't have both true
@@ -72,13 +75,10 @@
PrintStep := 5; % Coarseness, in resolvec
- defaultscale := 0.75;
- defaultfont := "cmss17"; % This is used by cartaxes
-
PageHeight := 9in;
PageWidth := 6in; % And this is used by produce_auto_scale
- MaxFearLimit := 15; % Valid Maximum Distance from Origin
+ MaxFearLimit := 17; % Valid Maximum Distance from Origin
HigColor := 0.85white; % These two colors are used in
SubColor := 0.35white; % fillfacewithlight
@@ -115,6 +115,8 @@
NCL := 0; % closedline
ForePen := pencircle scaled 15pt;
BackPen := pencircle scaled 9pt;
+
+ MalcomX := false;
%%% The variables PhotoMarks, PhotoPair[], PhotoPoint[]
%%% and CLPath[] have NO default values.
@@ -124,26 +126,42 @@
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%% Very basic:
-% Colors have three coordinates. Get one.
+% Colors have three or four coordinates. Get one.
- def X(expr A) =
- redpart A
+ def X(expr A) =
+ if color A: redpart A elseif MalcomX: blackpart A else: cyanpart A fi
enddef;
-
- def Y(expr A) =
- greenpart A
+
+ def Y(expr A) =
+ if color A: greenpart A else: magentapart A fi
enddef;
-
+
def Z(expr A) =
- bluepart A
+ if color A: bluepart A else: yellowpart A fi
enddef;
-
+
+ def W(expr A) =
+ blackpart A
+ enddef;
+
% The length of a vector.
def conorm(expr A) =
- ( X(A) ++ Y(A) ++ Z(A) )
+ ( X(A) ++ Y(A) ++ Z(A) )
enddef;
+ def cmyknorm(expr A) = %% This is not good when MalcomX is true
+ ( X(A) ++ Y(A) ++ Z(A) ++ W(A) )
+ enddef;
+
+ def makecmyk( expr A, B ) =
+ ( ( X(A), Y(A), Z(A), B ) )
+ enddef;
+
+ def maketrio( expr A ) =
+ ( ( X(A), Y(A), Z(A) ) )
+ enddef;
+
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%% Vector Calculus:
@@ -222,6 +240,13 @@
endgroup
enddef;
+% The following routine is used by circularsheet and may be used to
+% rotate vectors elliptically.
+
+ vardef planarrotation( expr VecX, VecY, TheAngle ) =
+ ( VecX*cosd( TheAngle ) + VecY*sind( TheAngle ) )
+ enddef;
+
% The following routine could be used by kindofcube and may be used to
% rotate polyhedra (must cycle through all Vs before calling makeface).
@@ -240,6 +265,20 @@
endgroup
enddef;
+% Rotate a vector around another. Supposes all vectors share the same origin.
+
+ def rotvecaroundanother( expr Angle, RotVec, FixVec ) =
+ begingroup
+ save uf, cf, xr, yr;
+ color uf, cf, xr, yr;
+ uf = N( FixVec );
+ yr = ccrossprod( uf, RotVec );
+ cf = uf*cdotprod( uf, RotVec );
+ xr = RotVec - cf;
+ ( cf + planarrotation( xr, yr, Angle ) )
+ endgroup
+ enddef;
+
% inplanarvolume is used by kindofcube.
def inplanarvolume( expr PointPerpA, PointPerpB, Point ) =
@@ -308,6 +347,10 @@
endgroup
enddef;
+% The following macro fits a figure in the page.
+% Probably it is obsolete since MetaPost 1.000
+% Should be the last command in a figure.
+
def produce_auto_scale =
begingroup
picture storeall, scaleall;
@@ -325,30 +368,29 @@
endgroup
enddef;
- vardef cstr( expr Cl ) =
- "(" &
- decimal(X(Cl)) &
- "," &
- decimal(Y(Cl)) &
- "," &
- decimal(Z(Cl)) &
- ")"
- enddef;
-
- vardef bstr( expr bv ) =
- save bstring; string bstring;
- if bv: bstring = "true"; else: bstring = "false"; fi;
- bstring
- enddef;
+ vardef cstr( expr Cl ) =
+ "(" &
+ decimal(X(Cl)) &
+ "," &
+ decimal(Y(Cl)) &
+ "," &
+ decimal(Z(Cl)) &
+ ")"
+ enddef;
+
+ vardef bstr( expr bv ) =
+ save bstring; string bstring;
+ if bv: bstring = "true"; else: bstring = "false"; fi;
+ bstring
+ enddef;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%% Fundamental:
% Rigorous Projection. This the kernel of all these lines of code.
% It won't work if R belongs the plane that contains f and that is
-% ortogonal to vector f, unless SphericalDistortion is true.
-% f must not be on a line parallel to zz and that contains the
-% viewcentr.
+% ortogonal to vector f-viewcentr, unless SphericalDistortion is true.
+% f-viewcentr must not be (anti-)parallel to zz.
def rp(expr R) =
begingroup
@@ -383,8 +425,8 @@
squarf = cdotprod( f-viewcentr, f-viewcentr );
radio = cdotprod( R-viewcentr, f-viewcentr );
eta = 1 - radio/squarf;
- if eta < 0.03: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% DANGER
- eta := 0.03;
+ if abs((horiz,verti)) > MaxFearLimit*eta:
+ eta := abs((horiz,verti))/MaxFearLimit;
fi;
fi;
projpoi = (horiz,verti)/eta;
@@ -467,7 +509,7 @@
endgroup
enddef;
- def drawsegment( expr A, B )=
+ def drawsegment( expr A, B ) =
begingroup
if SphericalDistortion:
draw pathofstraightline( A, B );
@@ -496,6 +538,24 @@
endgroup
enddef;
+% Orthogonal axes with prescribed lengths and labels.
+
+ def orthaxes(expr axex, strx, axey, stry, axez, strz ) =
+ begingroup
+ save axxc, ayyc, azzc;
+ color axxc, ayyc, azzc;
+ axxc = (axex,0,0);
+ ayyc = (0,axey,0);
+ azzc = (0,0,axez);
+ drawarrow rp(black)..rp(axxc);
+ drawarrow rp(black)..rp(ayyc);
+ drawarrow rp(black)..rp(azzc);
+ label.bot( strx ,rp(axxc));
+ label.bot( stry ,rp(ayyc));
+ label.lft( strz ,rp(azzc));
+ endgroup
+ enddef;
+
% This is it. Draw an arch beetween two straight lines with a
% common point (Or) in three-dimensional-euclidian-space and
% place a label near the middle of the arch. Points A and B
@@ -579,9 +639,9 @@
def rigorouscircle( expr CenterPos, AngulMom, Radius ) =
begingroup
- save ind, G, Dna, Dnb, al;
+ save ind, G, Dn, Dna, Dnb, al, vec;
numeric ind, G;
- color vec[], Dna, Dnb;
+ color vec[], Dn, Dna, Dnb;
path al;
vec1 = ncrossprod( CenterPos-f, AngulMom);
for ind=2 step 2 until 8:
@@ -592,8 +652,9 @@
al = rp(Radius*vec1+CenterPos)
for ind=2 upto 8:
hide(
- Dna:=ncrossprod(ncrossprod(vec[ind-1],vec[ind]),vec[ind-1]);
- Dnb:=ncrossprod(ncrossprod(vec[ind],vec[ind-1]),vec[ind])
+ Dn:=ncrossprod(vec[ind-1],vec[ind]);
+ Dna:=ncrossprod(Dn,vec[ind-1]);
+ Dnb:=ncrossprod(-Dn,vec[ind])
)
..controls rp(Radius*vec[ind-1]+CenterPos+G*Dna)
and rp(Radius*vec[ind] +CenterPos+G*Dnb)
@@ -637,8 +698,7 @@
% Draw lines with a better expression of three-dimensionality.
- def emptyline(expr JoinP,ThickenFactor,OutCol,InCol,theN,EmptyFrac,sN)
- (text LinFunc) =
+ def emptyline(expr JoinP,ThickenFactor,OutCol,InCol,theN,EmptyFrac,sN)(text LinFunc) =
begingroup
save i, j, k;
numeric i, j, k;
@@ -667,8 +727,7 @@
% Draw space-paths of possibly closed lines making use of "getready"
- def closedline( expr ThisIsClosed, theN, ForeFrac, BackFrac )
- ( text LinFunc ) =
+ def closedline( expr ThisIsClosed, theN, ForeFrac, BackFrac )( text LinFunc ) =
begingroup
save i, comm;
numeric i;
@@ -748,8 +807,10 @@
% In order to complete the drawing of this solid you have to choose one of
% the edges to be drawn immediatly afterwards.
- def twocyclestogether( expr CycleA, CycleB )=
+ def twocyclestogether( expr CycleA, CycleB ) =
begingroup
+ save TheLengthOfA, TheLengthOfB, TheMargin, Leng, i;
+ save SubPathA, SubPathB, PolygonPath, FinalPath;
numeric TheLengthOfA, TheLengthOfB, TheMargin, Leng, i;
path SubPathA, SubPathB, PolygonPath, FinalPath;
TheMargin = 0.02;
@@ -772,19 +833,11 @@
def ellipticpath(expr CenterPos, OneAxe, OtherAxe ) =
begingroup
- save cirath, ind, vec;
+ save ind;
numeric ind;
- color vec[];
- path cirath;
- for ind=1 upto 36:
- vec[ind] = CenterPos+OneAxe*cosd(ind*10)+OtherAxe*sind(ind*10);
- endfor;
- cirath = rp( vec1 )
- for ind=2 upto 36:
- ...rp( vec[ind] )
- endfor
- ...cycle;
- ( cirath )
+ ( for ind=1 upto 36:
+ rp( CenterPos+planarrotation(OneAxe,OtherAxe,ind*10) )...
+ endfor cycle )
endgroup
enddef;
@@ -792,27 +845,18 @@
def ellipticshadowpath(expr CenterPos, OneAxe, OtherAxe ) =
begingroup
- save cirath, ind, vec;
+ save ind;
numeric ind;
- color vec[];
- path cirath;
- for ind=1 upto 36:
- vec[ind] = CenterPos+OneAxe*cosd(ind*10)+OtherAxe*sind(ind*10);
- endfor;
- cirath = rp( cb( vec1 ) )
- for ind=2 upto 36:
- ...rp( cb( vec[ind] ) )
- endfor
- ...cycle;
- ( cirath )
+ ( for ind=1 upto 36:
+ rp( cb( CenterPos+planarrotation(OneAxe,OtherAxe,ind*10) ) )...
+ endfor cycle )
endgroup
enddef;
% It should be possible to attach some text to some plan.
% Unfortunately, this only works correctly when ParallelProj := true;
- def labelinspace(expr KeepRatio,RefPoi,BaseVec,UpVec)
- (text SomeString)=
+ def labelinspace(expr KeepRatio,RefPoi,BaseVec,UpVec)(text SomeString) =
begingroup
save labelpic, plak, lrc, ulc, llc, centerc, aratio, newbase;
picture labelpic;
@@ -840,12 +884,13 @@
% It should be possible to attach some path to some surface.
- def closedpathinspace( expr SomeTDPath )( text ConverterFunc )=
+ def closedpathinspace( expr SomeTDPath, NDivide )( text ConverterFunc ) =
begingroup
- save i, outpath;
- numeric i;
+ save i, outpath, st;
+ numeric i, st;
path outpath;
- outpath = for i=0.25 step 0.25 until (length SomeTDPath):
+ st = 1/NDivide;
+ outpath = for i=st step st until (length SomeTDPath):
ConverterFunc( point i of SomeTDPath ) --
endfor cycle;
( outpath )
@@ -860,9 +905,9 @@
def goodcirclepath(expr CenterPos, AngulMom, Radius ) =
begingroup
- save cirath, vecx, vecy, ind, goodangulmom, decision;
+ save cirath, vecx, vecy, ind, goodangulmom, decision, view;
numeric ind, decision;
- color vecx, vecy, vec[], goodangulmom, view;
+ color vecx, vecy, goodangulmom, view;
path cirath;
view = f-CenterPos;
decision = cdotprod( view, AngulMom );
@@ -871,19 +916,13 @@
else:
goodangulmom = AngulMom;
fi;
- vecx = ncrossprod( view, goodangulmom );
+ vecx = Radius*ncrossprod( view, goodangulmom );
decision := getangle( view, goodangulmom );
if decision > 0.5: %%%%%%%%%%%%%%% DANGER %%%
- vecy = ncrossprod( goodangulmom, vecx );
- for ind=1 upto 36:
- vec[ind] := vecx*cosd(ind*10) + vecy*sind(ind*10);
- vec[ind] := CenterPos + vec[ind]*Radius;
- endfor;
- cirath = rp( vec1 )
- for ind=2 upto 36:
- ...rp( vec[ind] )
- endfor
- ...cycle;
+ vecy = Radius*ncrossprod( goodangulmom, vecx );
+ cirath = for ind=10 step 10 until 360:
+ rp( CenterPos + planarrotation(vecx,vecy,ind) )...
+ endfor cycle;
else:
cirath = head_on_circle( CenterPos, Radius );
fi;
@@ -919,23 +958,15 @@
def head_on_circle(expr Pos, Radius ) =
begingroup
- save cirath, vecx, vecy, ind, view;
+ save vecx, vecy, ind, view;
numeric ind;
- color vecx, vecy, vec[], view;
- path cirath;
+ color vecx, vecy, view;
view = f-Pos;
- vecx = N( (-Y(view), X(view), 0) );
- vecy = ncrossprod( view, vecx );
- for ind=1 upto 36:
- vec[ind] := vecx*cosd(ind*10) + vecy*sind(ind*10);
- vec[ind] := Pos + vec[ind]*Radius;
- endfor;
- cirath = rp( vec1 )
- for ind=2 upto 36:
- ...rp( vec[ind] )
- endfor
- ...cycle;
- ( cirath )
+ vecx = Radius*N( (-Y(view), X(view), 0) );
+ vecy = Radius*ncrossprod( view, vecx );
+ ( for ind=10 step 10 until 360:
+ rp( Pos + planarrotation(vecx,vecy,ind) )...
+ endfor cycle )
endgroup
enddef;
@@ -945,7 +976,10 @@
def spatialhalfcircle(expr Center, AngulMom, Radius, ItsTheNearest ) =
begingroup
- save auxil, auxih;
+ save va, vb, vc, cc, vd, ux, uy, pa, pb;
+ save nr, cn, valx, valy, valr, choiceang;
+ save auxil, auxih, fcirc, returnp;
+ save choice;
color va, vb, vc, cc, vd, ux, uy, pa, pb;
numeric nr, cn, valx, valy, valr, choiceang;
path auxil, auxih, fcirc, returnp;
@@ -957,7 +991,7 @@
vd := cc - Center; % vd := va + vc;
nr := conorm( vd );
if Radius >= nr:
- returnp := rp( cc );
+ returnp := rp( cc ); % this single point will show up in spheroid
else:
valr := Radius*Radius;
valx := valr/nr;
@@ -1011,7 +1045,7 @@
begingroup
save va, vb, vc, cc, vd, base, holepic;
save vA, cC, nr, vala, valb, hashole, istube;
- save auxil, auxih, rect, halfl, halfh, thehole;
+ save auxil, auxih, halfl, halfh, thehole;
save auxili, auxihi, rect, theshadow;
color va, vb, vc, cc, vd, base;
@@ -1126,14 +1160,14 @@
% draw only the in fact visible part of circular lines. Please, don't
% put the vertex too close to the base plan when UsualForm=false.
- def rigorouscone(expr UsualForm,CenterPos,AngulMom,Radius,VertexPos)=
+ def rigorouscone(expr UsualForm,CenterPos,AngulMom,Radius,VertexPos) =
begingroup
- save basepath, themargin, thelengthofc, thesubpath, fullpath;
- save newlen, finalpath, i, auxpath, pa, pb, auxt, bigcirc;
- save startt, endt;
+ save basepath, thesubpath, fullpath, finalpath, auxpath, bigcirc;
+ save themargin, newlen, i, auxt, startt, endt, thelengthofc;
+ save pa, pb, pc, pd, pe;
path basepath, thesubpath, fullpath, finalpath, auxpath;
path bigcirc;
- numeric themargin, newlen, i, auxt, startt, endt;
+ numeric themargin, newlen, i, auxt, startt, endt, thelengthofc;
pair pa, pb, pc, pd, pe;
themargin = 0.02; %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% DANGER
basepath = goodcirclepath( CenterPos, AngulMom, Radius );
@@ -1169,7 +1203,7 @@
endgroup
enddef;
- def verygoodcone(expr BackDash,CenterPos,AngulMom,Radius,VertexPos)=
+ def verygoodcone(expr BackDash,CenterPos,AngulMom,Radius,VertexPos) =
begingroup
save thepath, lenpath, bonevec, sidevec, viewaxe, cipath;
save thelengthofc, thesubpath, themargin, basepath;
@@ -1219,7 +1253,7 @@
def rigorousfearpath(expr Center, Radius ) =
begingroup
- save auxil;
+ save auxil, ux, uy, newcen, nr, valx, valy, valr;
color ux, uy, newcen;
numeric nr, valx, valy, valr;
path auxil;
@@ -1235,7 +1269,10 @@
def rigorousfearshadowpath(expr Center, Radius ) =
begingroup
- save auxil, auxih;
+ save ux, uy, newcen;
+ save nr, valx, valy, valr, lenr;
+ save auxil, auxih, fcirc, returnp;
+ save dcenter;
color ux, uy, newcen;
numeric nr, valx, valy, valr, lenr;
path auxil, auxih, fcirc, returnp;
@@ -1252,7 +1289,7 @@
% It's a globe (without land).
- def tropicalglobe( expr NumLats, TheCenter, Radius, AngulMom )=
+ def tropicalglobe( expr NumLats, TheCenter, Radius, AngulMom ) =
begingroup
save viewaxe, sinalfa, sinbeta, globaxe, aux, limicos, lc;
save stepang, actang, newradius, foc, newcenter, cpath, i;
@@ -1315,7 +1352,7 @@
% An elliptical frustum:
- def whatisthis(expr CenterPos,OneAxe,OtherAxe,CentersDist,TheFactor)=
+ def whatisthis(expr CenterPos,OneAxe,OtherAxe,CentersDist,TheFactor) =
begingroup
save patha, pathb, pathc, centersvec, noption;
path patha, pathb, pathc;
@@ -1349,17 +1386,122 @@
fi
endgroup
enddef;
-
+
+% Probably the last algorithm I'm going to write for featpost...
+
+ def spheroidshadow( expr CentrPoi, NorthPoleVec, Ray ) =
+ begingroup
+ save a, k, fx, fy, tmpa, tmpb, tmpc, ep, vax, wax, xax, yax, zax, cm, cp;
+ save bdh, bdm, bdp, bdv, sm, sp, i, cac;
+ numeric a, k, fx, fy, tmpa, tmpb, tmpc, cm, cp, sm, sp, i;
+ path ep;
+ color vax, wax, xax, yax, zax, cac;
+ pair bdh, bdm, bdp, bdv;
+ vax = LightSource-CentrPoi;
+ if cdotprod(NorthPoleVec,vax)<0:
+ xax = -N(NorthPoleVec);
+ else:
+ xax = N(NorthPoleVec);
+ fi;
+ a = conorm(NorthPoleVec);
+ k = a/Ray;
+ if getangle(xax,vax) > 0.5: %%%%%%%%%%%%%%% DANGER %%%
+ zax = ncrossprod(xax,vax);
+ else:
+ zax = N( ( 0, Z(vax), -Y(vax) ) );
+ fi;
+ yax = ncrossprod(zax,xax);
+ fx = cdotprod(vax,xax);
+ fy = cdotprod(vax,yax);
+ tmpa = Ray*fx/k;
+ tmpb = fy*(fy ++ (fx/k) +-+ Ray);
+ tmpc = ((fx/k)**2)+(fy**2);
+ cm = (tmpa-tmpb)/tmpc;
+ cp = (tmpa+tmpb)/tmpc;
+ sm = 1 +-+ cm;
+ if fx<a:
+ sp = 1 +-+ cp;
+ else:
+ sp = -( 1 +-+ cp );
+ fi;
+ bdm = (k*cm,sm)*Ray;
+ bdp = (k*cp,sp)*Ray;
+ bdh = 0.5[bdp,bdm];
+ tmpc := Ray*( 1 +-+ ((xpart bdh)/a) );
+ tmpb := tmpc +-+ (ypart bdh);
+ bdv = bdm-bdp;
+ wax = 0.5*( xax*( xpart (bdv) ) + yax*( ypart (bdv) ) );
+ cac = CentrPoi+ xax*( xpart (bdh) ) + yax*( ypart (bdh) );
+ fill ellipticshadowpath( cac, wax, zax*tmpb );
+ endgroup
+ enddef;
+
+ def spheroid( expr CentrPoi, NorthPoleVec, Ray ) =
+ begingroup
+ save a, k, fx, fy, tmpa, tmpb, tmpc, ep;
+ save vax, wax, xax, yax, zax, cm, cp;
+ save bdh, bdm, bdp, bdv, sm, sp, i, cac;
+ numeric a, k, fx, fy, tmpa, tmpb, tmpc, cm, cp, sm, sp, i;
+ path ep;
+ color vax, wax, xax, yax, zax, cac;
+ pair bdh, bdm, bdp, bdv;
+ if ShadowOn:
+ spheroidshadow( CentrPoi, NorthPoleVec, Ray );
+ fi;
+ vax = f-CentrPoi;
+ if cdotprod(NorthPoleVec,vax)<0:
+ xax = -N(NorthPoleVec);
+ else:
+ xax = N(NorthPoleVec);
+ fi;
+ a = conorm(NorthPoleVec);
+ k = a/Ray;
+ if getangle(xax,vax) > 0.5: %%%%%%%%%%%%%%% DANGER %%%
+ zax = ncrossprod(xax,vax);
+ else:
+ zax = N( (-Y(vax), X(vax), 0) );
+ fi;
+ yax = ncrossprod(zax,xax);
+ fx = cdotprod(vax,xax);
+ fy = cdotprod(vax,yax);
+ tmpa = Ray*fx/k;
+ tmpb = fy*(fy ++ (fx/k) +-+ Ray);
+ tmpc = ((fx/k)**2)+(fy**2);
+ cm = (tmpa-tmpb)/tmpc;
+ cp = (tmpa+tmpb)/tmpc;
+ sm = 1 +-+ cm;
+ if fx<a:
+ sp = 1 +-+ cp;
+ else:
+ sp = -( 1 +-+ cp );
+ fi;
+ bdm = (k*cm,sm)*Ray;
+ bdp = (k*cp,sp)*Ray;
+ bdh = 0.5[bdp,bdm];
+ tmpc := Ray*( 1 +-+ ((xpart bdh)/a) );
+ tmpb := tmpc +-+ (ypart bdh);
+ bdv = bdm-bdp;
+ wax = 0.5*( xax*( xpart (bdv) ) + yax*( ypart (bdv) ) );
+ cac = CentrPoi+ xax*( xpart (bdh) ) + yax*( ypart (bdh) );
+ ep = ellipticpath( cac, wax, zax*tmpb );
+ unfill ep;
+ draw ep;
+ draw spatialhalfcircle( CentrPoi, NorthPoleVec, Ray, true );
+ 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.
- def fakehole( expr CenterPos, LenVec, Radius )=
+ def fakehole( expr CenterPos, LenVec, Radius ) =
begingroup
save patha, pathb, pathc, noption, hashole, auxv, poption, vv;
- path patha, pathb, pathc;
- numeric noption;
- boolean hashole, poption;
+ save ta, tb, taf, tbf, margint, stopair, pa, pb, testpath, isin;
+ path patha, pathb, pathc, pa, pb, testpath;
+ numeric noption, ta, tb, margint;
+ boolean hashole, poption, isin;
color auxv, vv;
+ pair stopair;
vv = f-CenterPos;
patha := rigorouscircle( CenterPos, LenVec, Radius );
pathb := rigorouscircle( CenterPos+LenVec, LenVec, Radius );
@@ -1371,22 +1513,35 @@
draw patha;
draw pathb;
else:
-% draw patha withcolor green; show patha;
-% draw pathb withcolor green; show pathb;
- hashole := (-1,-1) <> ( patha intersectiontimes pathb );
- if hashole:
- pathc := buildcycle( patha, pathb ); % I don't get it!
- %fill pathc withcolor red; % see fakehole.mp
- fi;
noption = cdotprod( LenVec, vv );
if noption > (conorm(LenVec)**2):
- draw pathb;
- if hashole:
- draw pathc;
- fi;
+ pa = patha;
+ pb = pathb;
elseif noption < 0:
- draw patha;
- if hashole:
+ pa = pathb;
+ pb = patha;
+ fi;
+ draw pb;
+ stopair = pa intersectiontimes pb;
+ hashole = (-1,-1) <> stopair;
+ if hashole:
+ testpath = rp(CenterPos+0.5*LenVec)--(point 0 of pa);
+ isin = (-1,-1) <> testpath intersectiontimes pb;
+ if not isin:
+ ta = xpart stopair;
+ tb = ypart stopair;
+ stopair := (reverse pa) intersectiontimes (reverse pb);
+ taf = length pa - xpart stopair;
+ tbf = length pb - ypart stopair;
+ margint = 0.01; % DANGER!
+ draw (subpath (0,ta-margint) of pa)--
+ (subpath (tb+margint,tbf-margint) of pb)--
+ (subpath (taf+margint,length pa - margint) of pa)--
+ cycle;
+ else:
+ pathc := buildcycle( pa, pb ); % I don't get it!
+ % Why doesn't buildcycle work all the time??? See fakehole.mp
+ % When point 0 of pa is inside pb, builcycle doesn't work!!
draw pathc;
fi;
fi;
@@ -1396,8 +1551,7 @@
% It is time for a kind of cube. Don't use SphericalDistortion here.
- def kindofcube(expr WithDash, IsVertex, RefP,
- AngA, AngB, AngC, LenA, LenB, LenC ) =
+ def kindofcube(expr WithDash, IsVertex, RefP, AngA, AngB, AngC, LenA, LenB, LenC ) =
begingroup
save star, pos, patw, patb, refv, near, centre, farv;
save newa, newb, newc, veca, vecb, vecc, auxx, auxy, i;
@@ -1498,7 +1652,7 @@
% It's a bit late now but the stage must be set.
- def setthestage( expr NumberOfSideSquares, SideSize )=
+ def setthestage( expr NumberOfSideSquares, SideSize ) =
begingroup
save i, j, squaresize, squarepath, ca, cb, cc, cd;
numeric i, j, squaresize;
@@ -1519,7 +1673,7 @@
endgroup
enddef;
- def setthearena( expr NumberOfDiameterCircles, ArenaDiameter )=
+ def setthearena( expr NumberOfDiameterCircles, ArenaDiameter ) =
begingroup
save i, j, circlesize, polar, currpos, phi, cpath;
numeric i, j, circlesize, polar, phi;
@@ -1718,6 +1872,88 @@
endgroup
enddef;
+% Take a "quarter" of a "stretched" donut (under construction)
+
+ def quartertorus( expr Tcenter, Tstart, Tfinis, Sray ) =
+ begingroup
+ save sideaxe, viewline, circlecenter, circlemoment;
+ save i, angstep, cuspcond, coofrac, tmoment, tstart, tfinis;
+ save cpath, opath, ipath, wp, ep, refpair, fpath;
+ save vstart, vfinis, ostart, ofinis, cstart, cfinis;
+ color sideaxe, viewline, circlecenter, circlemoment;
+ color tmoment, vstart, vfinis, tstart, tfinis;
+ numeric i, angstep, coofrac;
+ path cpath, opath, ipath, wp, ep, cstart, cfinis, fpath;
+ pair outerp[], innerp[], refpair;
+ boolean cuspcond, ostart, ofinis;
+ angstep = 6;
+ viewline = f-Tcenter;
+% if cdotprod( viewline, Tstart ) < cdotprod( viewline, Tfinis ):
+% tstart = Tfinis;
+% tfinis = Tstart;
+% else:
+ tstart = Tstart;
+ tfinis = Tfinis;
+% fi;
+ tmoment = ncrossprod( tstart, tfinis );
+ vstart = ncrossprod( tstart, tmoment );
+ vfinis = ncrossprod( tmoment, tfinis );
+ ostart = cdotprod( viewline-tstart, vstart ) > 0;
+ ofinis = cdotprod( viewline-tfinis, vfinis ) > 0;
+ cstart = spatialhalfcircle(tstart+Tcenter,vstart,Sray,true);
+ cfinis = spatialhalfcircle(tfinis+Tcenter,vfinis,Sray,true);
+ if cdotprod( viewline, tmoment ) < 0:
+ tmoment := -tmoment;
+ fi;
+ refpair = unitvector( rp(Tcenter+tmoment)-rp(Tcenter) );
+ sideaxe = ncrossprod( tmoment, viewline );
+ coofrac = cdotprod( viewline, tmoment )/Sray;
+
+ for i=0 step angstep until 90:
+ circlecenter:= tstart*cosd(i)+tfinis*sind(i);
+ circlemoment:= ccrossprod(circlecenter,tmoment);
+ cpath:=spatialhalfcircle(circlecenter+Tcenter,circlemoment,Sray,true);
+ if cdotprod( sideaxe, circlecenter ) < 0:
+ outerp[i/angstep]=point 0 of cpath;
+ innerp[i/angstep]=point (length cpath) of cpath;
+ else:
+ innerp[i/angstep]=point 0 of cpath;
+ outerp[i/angstep]=point (length cpath) of cpath;
+ fi;
+ endfor;
+ opath = outerp0 for i=angstep step angstep until 90:
+ ..outerp[i/angstep] endfor;
+ ipath = innerp0 for i=angstep step angstep until 90:
+ ..innerp[i/angstep] endfor;
+ fpath = ipath---cfinis---reverse opath---cstart---cycle;
+ unfill fpath;
+ draw fpath;
+% draw cstart---ipath---cfinis;
+ if ostart:
+ cpath := rigorouscircle( Tcenter+tstart, vstart, Sray );
+ unfill cpath;
+ draw cpath;
+ fi;
+ if ofinis:
+ cpath := rigorouscircle( Tcenter+tfinis, vfinis, Sray );
+ unfill cpath;
+ draw cpath;
+ fi;
+
+% i := 0;
+% cuspcond = false;
+% forever:
+% exitif i > (90/angstep)-2;
+% i := incr( i );
+% cuspcond :=
+% refpair dotprod outerp[i+1] <
+% refpair dotprod outerp[i];
+% exitif cuspcond;
+% endfor;
+% undraw subpath (i,90/angstep) of opath;
+ endgroup
+ enddef;
+
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%% Non-standard objects:
@@ -1747,7 +1983,7 @@
endgroup
enddef;
- def simplecar(expr RefP, AngCol, LenCol, FronWheelCol, RearWheelCol )=
+ def simplecar(expr RefP, AngCol, LenCol, FronWheelCol, RearWheelCol ) =
begingroup
save veca, vecb, vecc, anga, angb, angc, lena, lenb, lenc;
save auxn, viewline, auxm, fl, fr, rl, rr, auxx, auxy;
@@ -1835,8 +2071,9 @@
% Oh! Well... I couldn't do without differential equations.
% The point is that I want to draw vectorial field lines in space.
% Keep it simple: second-order Runge-Kutta method.
+% This is for solving first order differential equations
- def fieldlinestep( expr Spos, Step )( text VecFunc )=
+ def fieldlinestep( expr Spos, Step )( text VecFunc ) =
begingroup
save kone, ktwo;
color kone, ktwo;
@@ -1846,7 +2083,7 @@
endgroup
enddef;
- def fieldlinepath( expr Numb, Spos, Step )( text VecFunc )=
+ def fieldlinepath( expr Numb, Spos, Step )( text VecFunc ) =
begingroup
save ind, flpath, prevpos, thispos;
numeric ind;
@@ -1864,8 +2101,9 @@
enddef;
% Another point is that I want to draw trajectories in space.
+% This is for solving second order differential equations
- def trajectorypath( expr Numb, Spos, Svel, Step )( text VecFunc )=
+ def trajectorypath( expr Numb, Spos, Svel, Step )( text VecFunc ) =
begingroup
save ind, flpath, prevpos, thispos, prevvel, thisvel;
save rone, rtwo, vone, vtwo;
@@ -1895,10 +2133,47 @@
endgroup
enddef;
+% Another point is that I want to draw trajectories in space and
+% dependant on velocity: VecFunc( position, velocity ).
+% This time is fourth-order Runge-Kutta.
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% CHANGES PrintStep!!!!
+
+ def dragtrajectorypath( expr Spos, Svel, Step )( text VecFunc ) =
+ begingroup
+ save ind, flpath, prevpos, thispos, prevvel, thisvel;
+ save rone, rtwo, rthr, rfou, vone, vtwo, vthr, vfou;
+ numeric ind;
+ color prevpos, thispos, prevvel, thisvel;
+ color rone, rtwo, rthr, rfou, vone, vtwo, vthr, vfou;
+ path flpath;
+ prevpos = Spos;
+ prevvel = Svel;
+ flpath = rp( Spos );
+ ind = 1;
+ forever:
+ vone := Step*VecFunc( prevpos , prevvel );
+ rone := Step*prevvel;
+ vtwo := Step*VecFunc( prevpos+0.5*rone, prevvel+0.5*vone );
+ rtwo := Step*( prevvel+0.5*vone );
+ vthr := Step*VecFunc( prevpos+0.5*rtwo, prevvel+0.5*vtwo );
+ rthr := Step*( prevvel+0.5*vtwo );
+ vfou := Step*VecFunc( prevpos+rthr, prevvel+vthr );
+ rfou := Step*( prevvel+vthr );
+ thisvel := prevvel+(vtwo+vthr)/3+(vone+vfou)/6;
+ thispos := prevpos+(rtwo+rthr)/3+(rone+rfou)/6;
+ exitif Z( thispos ) < -0.0001; %%%%%%%%%% EDIT!
+ prevpos := thispos;
+ prevvel := thisvel;
+ flpath := flpath--rp( thispos );
+ endfor;
+ PrintStep := Y(thispos);
+ ( flpath )
+ endgroup
+ enddef;
+
% And now i stop.
- def magnetictrajectorypath( expr Numb, Spos, Svel, Step )
- ( text VecFunc )=
+ def magnetictrajectorypath( expr Numb, Spos, Svel, Step )( text VecFunc ) =
begingroup
save ind, flpath, prevpos, thispos, prevvel, thisvel;
save rone, rtwo, rthr, rfou, vone, vtwo, vthr, vfou;
@@ -1940,36 +2215,36 @@
%%%%%%%%%%%%%%%%%%%%%%%%%%%% Advanced 3D-Object Definition Functions %%%%%
% Please check the examples in planpht.mp or the default object below %%%%
- vardef makeline@#( text vertices )=
- save counter;
- numeric counter;
- counter = 0;
- for ind=vertices:
- counter := incr( counter );
- L@#p[counter] := V[ind];
- endfor;
- npl@# := counter;
- NL := @#
- enddef;
-
- vardef makeface@#( text vertices )=
- save counter;
- numeric counter;
- counter = 0;
- for ind=vertices:
- counter := incr( counter );
- F@#p[counter] := V[ind];
- endfor;
- npf@# := counter;
- NF := @#;
- FCD[NF] := false
- enddef;
-
- vardef getready( expr commstr, refpoi ) =
- Nobjects := incr( Nobjects );
- ostr[Nobjects] := commstr;
- RefDist[Nobjects] := conorm( f - refpoi )
- enddef;
+ vardef makeline@#( text vertices ) =
+ save counter;
+ numeric counter;
+ counter = 0;
+ for ind=vertices:
+ counter := incr( counter );
+ L@#p[counter] := V[ind];
+ endfor;
+ npl@# := counter;
+ NL := @#
+ enddef;
+
+ vardef makeface@#( text vertices ) =
+ save counter;
+ numeric counter;
+ counter = 0;
+ for ind=vertices:
+ counter := incr( counter );
+ F@#p[counter] := V[ind];
+ endfor;
+ npf@# := counter;
+ NF := @#;
+ FCD[NF] := false
+ enddef;
+
+ vardef getready( expr commstr, refpoi ) =
+ Nobjects := incr( Nobjects );
+ ostr[Nobjects] := commstr;
+ RefDist[Nobjects] := conorm( f - refpoi )
+ enddef;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% Definition of a 3D-Object
% define vertices
@@ -2294,7 +2569,7 @@
% FillDraw a face
- def face_invisible( expr Facen )( text LineAtribs )=
+ def face_invisible( expr Facen )( text LineAtribs ) =
begingroup
save ghost;
path ghost;
@@ -2363,7 +2638,7 @@
color somepo;
numeric fract;
(b-a)*fract = somepo-a;
- cdotprod(swapc,somepo)=cdotprod(swapc,f);
+ cdotprod(swapc,somepo) =cdotprod(swapc,f);
if (fract>0) and (fract<1):
crosspoin[counter] := somepo;
counter := incr(counter);
@@ -2442,7 +2717,7 @@
% Draw only the faces, rigorously projecting the edges.
- def faceraytrace(expr Press, Col)=
+ def faceraytrace(expr Press, Col) =
begingroup
save i, j, a, b;
numeric i, j;
@@ -2463,7 +2738,7 @@
% Fast test for your three-dimensional object
- def draw_all_test( expr AlsoDrawLines )=
+ def draw_all_test( expr AlsoDrawLines ) =
begingroup
save i, j, a, b;
numeric i, j;
@@ -2507,7 +2782,7 @@
% Don't use SphericalDistortion here.
- def fill_faces( text LineAtribs )=
+ def fill_faces( text LineAtribs ) =
begingroup
save i;
numeric i;
@@ -2594,8 +2869,7 @@
endgroup
enddef;
- def generateonebiax(expr Lin, Phi, Theta, Long,
- SndDirAngl, Base, Currpos ) =
+ def generateonebiax(expr Lin, Phi, Theta, Long, SndDirAngl, Base, Currpos ) =
begingroup
save basevec, longvec, u, v;
color basevec, longvec, u, v;
@@ -2615,7 +2889,7 @@
endgroup
enddef;
- def director_invisible( expr SortEmAll, ThickenFactor, CyclicLines )=
+ def director_invisible( expr SortEmAll, ThickenFactor, CyclicLines ) =
begingroup
save i, j, k, farone, thisfar;
save outerr, innerr, direc, ounum;
@@ -2681,6 +2955,123 @@
endgroup
enddef;
+% Now two routines to draw "bananas", well, sort of...
+% Initially coded by Pedro J. Sebastião
+
+ def circularsheet( expr CenterP, Rad, VecX, VecY, StartA, FinisA, Width ) =
+ begingroup
+ save vecz, ind;
+ color vecz;
+ numeric ind;
+ vecz = ncrossprod( VecX, VecY );
+ ( rp( CenterP + Width*vecz + Rad*planarrotation( VecX, VecY, StartA ) )
+ for ind=StartA+1 upto FinisA:
+ --rp( CenterP + Width*vecz + Rad*planarrotation( VecX, VecY, ind ))
+ endfor
+ for ind=FinisA downto StartA:
+ --rp( CenterP + Rad*planarrotation( VecX, VecY, ind ) )
+ endfor
+ --cycle )
+ endgroup
+ enddef;
+
+ def banana( expr CenterPos, Radius, AngleColor, Wid, Amp ) =
+ begingroup
+ save ind, sinbeta, cosbeta, aux, delta, angfvbx, angpos, angneg, au;
+ save bx, by, bz, fv, beta, gamma, alfa;
+ save outpath, outpathb, outpathc, visneg, vispos;
+ numeric ind, sinbeta, cosbeta, aux, delta, angfvbx, angpos, angneg, au;
+ numeric beta, gamma, alfa;
+ color bx, by, bz, fv;
+ path outpath, outpathb, outpathc;
+ boolean visneg, vispos;
+
+ alfa = X(AngleColor);
+ beta = Y(AngleColor);
+ gamma= Z(AngleColor);
+ bx = eulerrotation( alfa, beta, gamma, red );
+ by = eulerrotation( alfa, beta, gamma, green );
+ bz = eulerrotation( alfa, beta, gamma, blue );
+
+ au = cdotprod( f-CenterPos, by );
+ if 0 > au:
+ by := -by;
+ fi;
+ fv = cdotprod( f, bx )*bx + cdotprod( f, by )*by;
+ aux = conorm( fv-CenterPos );
+ if aux > Radius:
+ cosbeta = Radius/aux;
+ sinbeta = ( aux +-+ Radius )/aux;
+ delta = angle( cosbeta, sinbeta );
+ angfvbx = getangle( fv - CenterPos, bx );
+ if 0 > cdotprod( fv - CenterPos, by ):
+ angfvbx := -angfvbx;
+ fi;
+ angpos = delta + angfvbx;
+ angneg = angfvbx - delta;
+ if ( angneg > -Amp ) and ( Amp > angneg ):
+ visneg = true;
+ else:
+ visneg = false;
+ fi;
+ if ( angpos > -Amp ) and ( Amp > angpos ):
+ vispos = true;
+ else:
+ vispos = false;
+ fi;
+ if visneg and not vispos:
+ outpath = circularsheet(CenterPos,Radius,bx,by,-Amp,angneg,Wid);
+ unfill outpath;
+ draw outpath;
+ outpathb = circularsheet(CenterPos,Radius,bx,by,angneg,Amp,Wid);
+ unfill outpathb;
+ draw outpathb
+ fi;
+ if vispos and not visneg:
+ outpath = circularsheet(CenterPos,Radius,bx,by,-Amp,angpos,Wid);
+ unfill outpath;
+ draw outpath;
+ outpathb = circularsheet(CenterPos,Radius,bx,by,angpos,Amp,Wid);
+ unfill outpathb;
+ draw outpathb
+ fi;
+ if (not vispos) and (not visneg):
+ outpath = circularsheet(CenterPos,Radius,bx,by,-Amp,Amp,Wid);
+ unfill outpath;
+ draw outpath;
+ fi;
+ if vispos and visneg:
+ if 0 > cdotprod( f-CenterPos, bx ):
+ outpath=circularsheet(CenterPos,Radius,bx,by,angneg,angpos,Wid);
+ unfill outpath;
+ draw outpath;
+ outpathb = circularsheet(CenterPos,Radius,bx,by,-Amp,angneg,Wid);
+ unfill outpathb;
+ draw outpathb;
+ outpathc = circularsheet(CenterPos,Radius,bx,by,angpos,Amp,Wid);
+ unfill outpathc;
+ draw outpathc;
+ else:
+ outpathb = circularsheet(CenterPos,Radius,bx,by,-Amp,angneg,Wid);
+ unfill outpathb;
+ draw outpathb;
+ outpathc = circularsheet(CenterPos,Radius,bx,by,angpos,Amp,Wid);
+ unfill outpathc;
+ draw outpathc;
+ outpath=circularsheet(CenterPos,Radius,bx,by,angneg,angpos,Wid);
+ unfill outpath;
+ draw outpath;
+ fi;
+ fi;
+ else:
+ outpath = circularsheet( CenterPos, Radius, bx, by, -Amp, Amp, Wid );
+ unfill outpath;
+ draw outpath;
+ fi;
+ draw rp(CenterPos+Radius*bx)--rp(CenterPos+Radius*bx+Wid*bz);
+ endgroup
+ enddef;
+
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%% Plotting:
@@ -2817,7 +3208,7 @@
% ColAtribone the color used for drawing;
% Colordensity depends on distance of the face from f
- def face_drawfill( expr Facen, dmin_, dmax_ ,ColAtrib, ColAtribone )=
+ def face_drawfill( expr Facen, dmin_, dmax_ ,ColAtrib, ColAtribone ) =
begingroup
save j, ptmp, colfac_, coltmp_;
path ghost;
@@ -2854,7 +3245,7 @@
% ColAtrib=color for filling faces
% ColAtribone=color for drawing edges
- def draw_invisible( expr Option, DoJS, ColAtrib, ColAtribone )=
+ def draw_invisible( expr Option, DoJS, ColAtrib, ColAtribone ) =
begingroup
save i, j, a, b, thisfar, ptmp, farone;
numeric i, j, farone[], dist[], thisfar, distmin_, distmax_;
@@ -2945,11 +3336,7 @@
% Define parametric surfaces with a triangular mesh... unless a
% quadrangular mesh can do a fine, rigorous job just as well.
- def partrimesh( expr nt,ns,
- lowt,higt,lows,higs,
- lowx,higx,lowy,higy,
- lowz,higz,
- facz)( text parSurFunc ) =
+ def partrimesh( expr nt,ns,lowt,higt,lows,higs,lowx,higx,lowy,higy,lowz,higz,facz)( text parSurFunc ) =
begingroup
save i, j, k, posx, posy, posz;
save counter, stept, steps, poss, post, tmpaux;
@@ -3224,13 +3611,70 @@
( trypoi )
endgroup;
enddef;
+
+ 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):
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% Verify if a path is cyclic (written by Scott Pakin)
+
+ def is_cyclic expr cpath =
+ (point 0 of cpath = point (length cpath) of cpath)
+ enddef;
+
% Produce the schematics of a spring.
- def springpath( expr begp, endp, piturnum, piturnproj, spgfrac )=
+ def springpath( expr begp, endp, piturnum, piturnproj, spgfrac ) =
begingroup
boolean leftside;
numeric counter, springwidth;
@@ -3265,8 +3709,7 @@
% Summarize a great length in a zig-zag frontier line
- def zigzagfrontier( expr begp, endp, nzigs,
- dev, zthick, tthick, fthick, excol, incol )=
+ def zigzagfrontier( expr begp, endp, nzigs, dev, zthick, tthick, fthick, excol, incol ) =
begingroup
interim linecap := squared;
interim linejoin := mitered;
@@ -3302,7 +3745,7 @@
% The name says it all.
- def randomcirc( expr radi, stddev, numpois )=
+ def randomcirc( expr radi, stddev, numpois ) =
begingroup
numeric i, astep;
path ranc;
@@ -3377,195 +3820,491 @@
endgroup
enddef;
- def radialcross( expr A, la, B, lb, GoUp) =
- begingroup
- numeric x, y, xa, xb, ya, yb, YM, YA, La, Lb;
- numeric AA, BB, CC, auxil, na, nb, norm;
- pair As, Bs, selectedpoint;
- na = abs(A);
- nb = abs(B);
- norm := 0;
- for t = na, nb, la, lb:
- if norm < t:
- norm := t;
- fi;
- endfor;
- xa = xpart A/norm;
- xb = xpart B/norm;
- ya = ypart A/norm;
- yb = ypart B/norm;
- La = la/norm;
- Lb = lb/norm;
- if abs( ya - yb ) < 0.005 :
- x := La**2 - Lb**2 + xb**2 - xa**2;
- x := 0.5*x/( xb - xa );
- auxil := sqrt( La**2 - (xa-x)**2 );
- As = ( x, ya + auxil );
- Bs = ( x, ya - auxil );
- else:
- YM := (xb-xa)/(ya-yb);
- YA := Lb**2 - La**2 + xa**2 - xb**2;
- YA := 0.5*( YA - (ya-yb)**2 )/(ya-yb);
- AA := 1 + YM**2;
- BB := 2*( YM*YA - xa );
- CC := xa**2 - La**2 + YA**2;
- CC := sqrt( BB**2 - 4*AA*CC );
- x := -0.5*( BB + CC )/AA;
- y := YA + ya + YM*x;
- Bs = ( x, y );
- x := -0.5*( BB - CC )/AA;
- y := YA + ya + YM*x;
- As = ( x, y );
- fi;
- if ypart As > ypart Bs:
- if GoUp:
- selectedpoint = As;
- else:
- selectedpoint = Bs;
- fi;
- elseif ypart As = ypart Bs:
- if xpart As > xpart Bs:
- if GoUp:
- selectedpoint = As;
- else:
- selectedpoint = Bs;
- fi;
- else:
- if GoUp:
- selectedpoint = Bs;
- else:
- selectedpoint = As;
- fi;
- fi;
- else:
- if GoUp:
- selectedpoint = Bs;
- else:
- selectedpoint = As;
- fi;
- fi;
- ( norm*selectedpoint )
- endgroup
- enddef;
-
- def ropethread( expr Index ) =
- begingroup
- save aux;
- numeric aux;
- if Index > RopeColors:
- aux = 0;
+ def radialcross( expr A, la, B, lb, GoUp) =
+ begingroup
+ numeric x, y, xa, xb, ya, yb, YM, YA, La, Lb;
+ numeric AA, BB, CC, auxil, na, nb, norm;
+ pair As, Bs, selectedpoint;
+ na = abs(A);
+ nb = abs(B);
+ norm := 0;
+ for t = na, nb, la, lb:
+ if norm < t:
+ norm := t;
+ fi;
+ endfor;
+ xa = xpart A/norm;
+ xb = xpart B/norm;
+ ya = ypart A/norm;
+ yb = ypart B/norm;
+ La = la/norm;
+ Lb = lb/norm;
+ if abs( ya - yb ) < 0.005 :
+ x := La**2 - Lb**2 + xb**2 - xa**2;
+ x := 0.5*x/( xb - xa );
+ auxil := La +-+ (xa-x);
+ As = ( x, ya + auxil );
+ Bs = ( x, ya - auxil );
+ else:
+ YM := (xb-xa)/(ya-yb);
+ YA := Lb**2 - La**2 + xa**2 - xb**2;
+ YA := 0.5*( YA - (ya-yb)**2 )/(ya-yb);
+ AA := 1 + YM**2;
+ BB := 2*( YM*YA - xa );
+ CC := xa**2 - La**2 + YA**2;
+ CC := sqrt( BB**2 - 4*AA*CC );
+ x := -0.5*( BB + CC )/AA;
+ y := YA + ya + YM*x;
+ Bs = ( x, y );
+ x := -0.5*( BB - CC )/AA;
+ y := YA + ya + YM*x;
+ As = ( x, y );
+ fi;
+ if ypart As > ypart Bs:
+ if GoUp:
+ selectedpoint = As;
+ else:
+ selectedpoint = Bs;
+ fi;
+ elseif ypart As = ypart Bs:
+ if xpart As > xpart Bs:
+ if GoUp:
+ selectedpoint = As;
else:
- aux = Index;
+ selectedpoint = Bs;
fi;
- ( aux )
- endgroup
- enddef;
-
- def ropepattern( expr BasePath, RopeWidth, Nturns ) =
- begingroup
- save indturns, nmoves, indthread, movelen, turnlen, totlen;
- numeric indturns, nmoves, indthread, movelen, turnlen, totlen;
- save lenpos, timar, steplen, indstep, startdownc, startupcol;
- numeric lenpos, timar, steplen, startdownc, indstep;
- save actuc, actdc, stepwidth;
- numeric actuc, actdc, stepwidth, startupcol;
- save p;
- pair p[];
- save actcolor;
- color actcolor;
- nmoves = 2*(RopeColors+1);
- totlen = arclength BasePath;
- turnlen = totlen/Nturns;
- movelen = turnlen/nmoves;
- steplen = movelen/2;
- startdownc = 0;
- startupcol = RopeColors;
- stepwidth = RopeWidth/RopeColors;
- for indturns=0 upto Nturns-1:
- for indmove=0 upto nmoves-1:
- for indstep=0 upto 3:
- lenpos :=
- indturns*turnlen+indmove*movelen+indstep*steplen;
- timar := arctime lenpos of BasePath;
- p[indstep] := direction timar of BasePath rotated 90;
- p[indstep] := unitvector( p[indstep] );
- p[indstep+4] := point timar of BasePath;
- endfor;
- actdc := startdownc;
- for indthread=0 upto RopeColors:
- p8 := p5-p1*(0.5*RopeWidth-(indthread-0.5)*stepwidth);
- p9 := p4-p0*(0.5*RopeWidth-indthread*stepwidth);
- p10:= p5-p1*(0.5*RopeWidth-(indthread+0.5)*stepwidth);
- p11:= p6-p2*(0.5*RopeWidth-indthread*stepwidth);
- actcolor := TableC[RopeColorSeq[actdc]];
- fill p8--p9--p10--p11--cycle withcolor actcolor;
- actdc := ropethread( incr( actdc ) );
- endfor;
- startdownc := ropethread( incr( startdownc ) );
- actuc := startupcol;
- p9 := p5+p1*0.5*(RopeWidth+stepwidth);
- p10:= p6+p2*0.5*RopeWidth;
- p11:= p7+p3*0.5*(RopeWidth+stepwidth);
- actcolor := TableC[RopeColorSeq[actuc]];
- fill p9--p10--p11--cycle withcolor actcolor;
- actuc := ropethread( incr( actuc ) );
- for indthread=0 upto RopeColors-1:
- p8 := p6+p2*(0.5*RopeWidth-indthread*stepwidth);
- p9 := p5+p1*(0.5*RopeWidth-(indthread+0.5)*stepwidth);
- p10:= p6+p2*(0.5*RopeWidth-(indthread+1)*stepwidth);
- p11:= p7+p3*(0.5*RopeWidth-(indthread+0.5)*stepwidth);
- actcolor := TableC[RopeColorSeq[actuc]];
- fill p8--p9--p10--p11--cycle withcolor actcolor;
- actuc := ropethread( incr( actuc ) );
- endfor;
- p8 := p6-p2*0.5*RopeWidth;
- p9 := p5-0.5*p1*(RopeWidth+stepwidth);
- p11:= p7-0.5*p3*(RopeWidth+stepwidth);
- actcolor := TableC[RopeColorSeq[actuc]];
- fill p8--p9--p11--cycle withcolor actcolor;
- startupcol := ropethread( incr( startupcol ) );
- endfor;
- endfor
- endgroup
- enddef;
+ else:
+ if GoUp:
+ selectedpoint = Bs;
+ else:
+ selectedpoint = As;
+ fi;
+ fi;
+ else:
+ if GoUp:
+ selectedpoint = Bs;
+ else:
+ selectedpoint = As;
+ fi;
+ fi;
+ ( norm*selectedpoint )
+ endgroup
+ enddef;
+
+ def ropethread( expr Index ) =
+ begingroup
+ save aux;
+ numeric aux;
+ if Index > RopeColors:
+ aux = 0;
+ else:
+ aux = Index;
+ fi;
+ ( aux )
+ endgroup
+ enddef;
+
+ def ropepattern( expr BasePath, RopeWidth, Nturns ) =
+ begingroup
+ save indturns, nmoves, indthread, movelen, turnlen, totlen;
+ numeric indturns, nmoves, indthread, movelen, turnlen, totlen;
+ save lenpos, timar, steplen, indstep, startdownc, startupcol;
+ numeric lenpos, timar, steplen, startdownc, indstep;
+ save actuc, actdc, stepwidth;
+ numeric actuc, actdc, stepwidth, startupcol;
+ save p;
+ pair p[];
+ save actcolor;
+ color actcolor;
+ nmoves = 2*(RopeColors+1);
+ totlen = arclength BasePath;
+ turnlen = totlen/Nturns;
+ movelen = turnlen/nmoves;
+ steplen = movelen/2;
+ startdownc = 0;
+ startupcol = RopeColors;
+ stepwidth = RopeWidth/RopeColors;
+ for indturns=0 upto Nturns-1:
+ for indmove=0 upto nmoves-1:
+ for indstep=0 upto 3:
+ lenpos :=
+ indturns*turnlen+indmove*movelen+indstep*steplen;
+ timar := arctime lenpos of BasePath;
+ p[indstep] := direction timar of BasePath rotated 90;
+ p[indstep] := unitvector( p[indstep] );
+ p[indstep+4] := point timar of BasePath;
+ endfor;
+ actdc := startdownc;
+ for indthread=0 upto RopeColors:
+ p8 := p5-p1*(0.5*RopeWidth-(indthread-0.5)*stepwidth);
+ p9 := p4-p0*(0.5*RopeWidth-indthread*stepwidth);
+ p10:= p5-p1*(0.5*RopeWidth-(indthread+0.5)*stepwidth);
+ p11:= p6-p2*(0.5*RopeWidth-indthread*stepwidth);
+ actcolor := TableC[RopeColorSeq[actdc]];
+ fill p8--p9--p10--p11--cycle withcolor actcolor;
+ actdc := ropethread( incr( actdc ) );
+ endfor;
+ startdownc := ropethread( incr( startdownc ) );
+ actuc := startupcol;
+ p9 := p5+p1*0.5*(RopeWidth+stepwidth);
+ p10:= p6+p2*0.5*RopeWidth;
+ p11:= p7+p3*0.5*(RopeWidth+stepwidth);
+ actcolor := TableC[RopeColorSeq[actuc]];
+ fill p9--p10--p11--cycle withcolor actcolor;
+ actuc := ropethread( incr( actuc ) );
+ for indthread=0 upto RopeColors-1:
+ p8 := p6+p2*(0.5*RopeWidth-indthread*stepwidth);
+ p9 := p5+p1*(0.5*RopeWidth-(indthread+0.5)*stepwidth);
+ p10:= p6+p2*(0.5*RopeWidth-(indthread+1)*stepwidth);
+ p11:= p7+p3*(0.5*RopeWidth-(indthread+0.5)*stepwidth);
+ actcolor := TableC[RopeColorSeq[actuc]];
+ fill p8--p9--p10--p11--cycle withcolor actcolor;
+ actuc := ropethread( incr( actuc ) );
+ endfor;
+ p8 := p6-p2*0.5*RopeWidth;
+ p9 := p5-0.5*p1*(RopeWidth+stepwidth);
+ p11:= p7-0.5*p3*(RopeWidth+stepwidth);
+ actcolor := TableC[RopeColorSeq[actuc]];
+ fill p8--p9--p11--cycle withcolor actcolor;
+ startupcol := ropethread( incr( startupcol ) );
+ endfor;
+ endfor
+ endgroup
+ enddef;
- def firsttangencypoint( expr Path, Point, ResolvN ) =
- begingroup
- save auxp, i, cutp, va, vb;
- path auxp;
- numeric i;
- pair cutp, va, vb;
- auxp =
- hide( va := unitvector( point 0 of Path - Point );
- vb := unitvector( direction 0 of Path ); )
- ( paircrossprod( va, vb ), 0 )
- for i=1/ResolvN step 1/ResolvN until length Path:
+ def firsttangencypoint( expr Path, Point, ResolvN ) =
+ begingroup
+ save auxp, i, cutp, va, vb;
+ path auxp;
+ numeric i;
+ pair cutp, va, vb;
+ auxp =
+ hide( va := unitvector( point 0 of Path - Point );
+ vb := unitvector( direction 0 of Path ); )
+ ( paircrossprod( va, vb ), 0 )
+ for i=1/ResolvN step 1/ResolvN until length Path:
hide( va := unitvector( point i of Path - Point );
- vb := unitvector( direction i of Path ); )
- ...( paircrossprod( va, vb ), i )
+ vb := unitvector( direction i of Path ); )
+ ...( paircrossprod( va, vb ), i )
+ endfor;
+ cutp = auxp intersectionpoint ( origin--( 0, length Path ) );
+ ( point ( ypart cutp ) of Path )
+ endgroup
+ 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.
+
+ def lasermachine( expr DefinedPath, Beam, CosLimit ) =
+ begingroup
+ save patlen, j;
+ save apoi, bpoi, cpoi, dpoi, epoi;
+ save apat, bpat, cpat, dpat, a, b, c;
+ save anew, bnew, cnew;
+ save aang, bang, cang, dang;
+ save newp, pairvector, cou, val;
+ numeric patlen, j;
+ pair apoi, bpoi, cpoi, dpoi, epoi;
+ pair apat, bpat, cpat, dpat, a, b, c;
+ pair anew, bnew, cnew, pairvector[];
+ numeric aang, bang, cang, dang, cou, val;
+ path newp;
+ patlen = length DefinedPath;
+ cou = 0;
+ for j=0 upto patlen-1:
+ apoi := precontrol j of DefinedPath;
+ bpoi := point j of DefinedPath;
+ cpoi := postcontrol j of DefinedPath;
+ dpoi := precontrol j+1 of DefinedPath;
+ epoi := point j+1 of DefinedPath;
+ apat := apoi-bpoi;
+ bpat := bpoi-cpoi;
+ cpat := cpoi-dpoi;
+ dpat := dpoi-epoi;
+ aang := angle( apat );
+ bang := angle( bpat );
+ cang := angle( cpat );
+ dang := angle( dpat );
+ val := cosd( 0.5*(aang-bang) );
+ bnew := cpoi+Beam*dir( 90+0.5*(bang+cang) )/cosd( 0.5*(bang-cang) );
+ cnew := dpoi+Beam*dir( 90+0.5*(cang+dang) )/cosd( 0.5*(cang-dang) );
+ if ( val > 0 ) and
+ ( val < CosLimit ) and
+ ( Beam*sind( 0.5*(aang-bang) ) < 0 ):
+
+ a := bpoi+Beam*dir(90+aang);
+ b := a-Beam*dir(aang);
+ anew := bpoi+Beam*dir(90+bang);
+ c := anew+Beam*dir(bang);
+ pairvector[3*cou] = a;
+ pairvector[3*cou+1] = b;
+ pairvector[3*cou+2] = c;
+ cou := incr(cou);
+ pairvector[3*cou] = anew;
+ pairvector[3*cou+1] = bnew;
+ pairvector[3*cou+2] = cnew;
+ cou := incr(cou);
+ else:
+ anew := bpoi+Beam*dir( 90+0.5*(aang+bang) )/val;
+ pairvector[3*cou] = anew;
+ pairvector[3*cou+1] = bnew;
+ pairvector[3*cou+2] = cnew;
+ cou := incr(cou);
+ fi;
endfor;
- cutp = auxp intersectionpoint ( origin--( 0, length Path ) );
- ( point ( ypart cutp ) of Path )
- endgroup
- enddef;
+ newp = for j=0 upto cou-1:
+ pairvector[3*j]..controls pairvector[3*j+1] and pairvector[3*j+2]..
+ endfor cycle;
+ ( newp )
+ endgroup
+ enddef;
+
+% Move the starting point of a cyclic path along that path
+ def startahead( expr DefinedPath, JumpTime ) =
+ begingroup
+ save patlen, j;
+ save apoi, bpoi, cpoi;
+ save newp;
+ numeric patlen, j;
+ pair apoi, bpoi, cpoi;
+ path newp;
+ patlen = length DefinedPath;
+ newp = for j=0 upto patlen-1:
+ hide(
+ apoi := point JumpTime+j of DefinedPath;
+ bpoi := postcontrol JumpTime+j of DefinedPath;
+ cpoi := precontrol JumpTime+j+1 of DefinedPath;
+ )
+ apoi..controls bpoi and cpoi..
+ endfor
+ cycle;
+ ( newp )
+ endgroup
+ enddef;
+
+% In order to use a "lasermachine" one needs a single full outline.
+% One may have two somewhat concentric cyclic paths intersecting in
+% several points.
+% The next routine may help but first the paths must be adapted with
+% "startahead" and/or "reverse" so that they both rotate in the same
+% direction and they start on consecutive "lobes" (hard to explain).
+% Now pay attention: given the direction of rotation (clockwise or
+% counter-clockwise) the SecondPath must start BEFORE the FirstPath.
+% And another problem: there must be at least four intersection points.
+% Very nasty routine. All because of finispoi...
+
+ def crossingline( expr FirstPath, SecondPath, TimeTolerance ) =
+ begingroup
+ save m;
+ save its, finispoi;
+ save increm, fo, ma, tmpp, mastarter;
+ numeric m;
+ pair its, finispoi;
+ path increm, fo, ma, tmpp, mastarter;
+ m = TimeTolerance;
+ fo = FirstPath;
+ ma = SecondPath;
+ its := fo intersectiontimes ma;
+ increm := subpath (m, (xpart its) - m ) of fo;
+ mastarter = subpath (m, (ypart its) - m ) of ma;
+ finispoi = reverse ma intersectionpoint reverse fo;
+
+ forever:
+
+ fo := subpath ( (xpart its)+m, length fo ) of fo;
+ ma := subpath ( (ypart its)+m, length ma ) of ma;
+ its := ma intersectiontimes fo;
+ tmpp := subpath ( m, (xpart its)-m ) of ma;
+ increm := increm...tmpp;
+
+ fo := subpath ( (ypart its)+m, length fo ) of fo;
+ ma := subpath ( (xpart its)+m, length ma ) of ma;
+ its := fo intersectiontimes ma;
+ tmpp := subpath ( m, (xpart its)-m ) of fo;
+ increm := increm...tmpp;
+
+ exitif abs( point (xpart its) of fo - finispoi ) < m;
+ endfor;
+
+ fo := subpath ( (xpart its)+m, length fo ) of fo;
+ ma := (subpath ( (ypart its)+m, (length ma)-m ) of ma)...mastarter;
+ its := ma intersectiontimes fo;
+ tmpp := subpath ( m, (xpart its)-m ) of ma;
+ increm := increm...tmpp;
+
+ tmpp := subpath ( (ypart its)+m, (length fo)-m ) of fo;
+ ( increm...tmpp...cycle )
+ endgroup
+ enddef;
+
% Calculate path areas (contributed by Boguslaw Jackowski
% to the metapost mailing list)
- vardef segmentarea( expr Ps ) =
- save xa, xb, xc, xd, ya, yb, yc, yd;
- ( xa, 20ya ) = point 0 of Ps;
- ( xb, 20yb ) = postcontrol 0 of Ps;
- ( xc, 20yc ) = precontrol 1 of Ps;
- ( xd, 20yd ) = point 1 of Ps;
- ( xb - xa )*( 10ya + 6yb + 3yc + yd )
- + ( xc - xb )*( 4ya + 6yb + 6yc + 4yd )
- + ( xd - xc )*( ya + 3yb + 6yc + 10yd )
- enddef;
+ vardef segmentarea( expr Ps ) =
+ save xa, xb, xc, xd, ya, yb, yc, yd;
+ ( xa, 20ya ) = point 0 of Ps;
+ ( xb, 20yb ) = postcontrol 0 of Ps;
+ ( xc, 20yc ) = precontrol 1 of Ps;
+ ( xd, 20yd ) = point 1 of Ps;
+ ( xb - xa )*( 10ya + 6yb + 3yc + yd )
+ + ( xc - xb )*( 4ya + 6yb + 6yc + 4yd )
+ + ( xd - xc )*( ya + 3yb + 6yc + 10yd )
+ enddef;
+
+ vardef cyclicpatharea( expr P ) = % result = area of the interior
+ segmentarea(subpath (0,1) of P)
+ for t=1 upto length(P)-1: + segmentarea(subpath (t,t+1) of P) endfor
+ enddef;
+
+% Mark bidimensional angles (contributed by Palle Jørgensen
+% to the metapost mailing list)
+
+ vardef archangle@#( expr _p, _q, _s, archwidth ) text _t =
+ begingroup;
+ save _a, _b, _w, _arch, _halfangle, _label_origin;
+ ( _a, _b ) = _p intersectiontimes _q;
+ pair _w;
+ _w = whatever[
+ point _a of _p +
+ archwidth * unitvector direction _a of _p,
+ point _a of _p +
+ archwidth * unitvector direction _a of _p +
+ (ypart.direction _a of _p, -xpart.direction _a of _p)
+ ];
+ _w = whatever[point _b of _q, point _b of _q + direction _b of _q];
+ path _arch;
+ _arch = point _a of _p +
+ archwidth * unitvector direction _a of _p{
+ (if direction _a of _p dotprod unitvector direction _b of _q > 0:
+ 1
+ else:
+ -1
+ fi) *
+ ( _w - (point _a of _p + archwidth * unitvector direction _a of _p) )
+ }..point _b of _q +
+ archwidth * unitvector direction _b of _q;
+ draw _arch _t;
+ path _halfangle;
+ _halfangle = point _a of _p - 2*archwidth*
+ unitvector( direction _a of _p + direction _b of _q )--point _a of _p +
+ 2*archwidth*unitvector( direction _a of _p + direction _b of _q );
+ pair _label_origin;
+ _label_origin = _halfangle intersectionpoint _arch;
+ label@#( _s, _label_origin ) _t;
+ endgroup;
+enddef;
- vardef cyclicpatharea( expr P ) = % result = area of the interior
- segmentarea(subpath (0,1) of P)
- for t=1 upto length(P)-1: + segmentarea(subpath (t,t+1) of P) endfor
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%%%% Part VI (strictly two-dimensional and related to planifications):
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% Default Values %%%%%%%%%%%%
+ pen grossa, usual, fina;
+ grossa= pencircle scaled 0.9mm;
+ usual = pencircle scaled 0.6mm;
+ fina = pencircle scaled 0.3mm;
+ numeric u, arestaref, defaultaresta;
+ u = 0.5mm;
+ defaultaresta = 30mm;
+ arestaref = defaultaresta;
+ picture pequeno, grande, nada;
+ linecap := squared; %%%%%%%%%%%%%%%%%%%%%
+ pequeno = dashpattern( off 3u on 1u );
+ grande = dashpattern( off 8u on 6u );
+ nada = dashpattern( on 1pt );
+ linecap := rounded; %%%%%%%%%%%%%%%%%%%%%
+ pair urCorner,ulCorner,llCorner,lrCorner,PageCenter;
+ urCorner := (8.13in,11.533in);
+ llCorner := (0.118in,35.5bp);
+ ulCorner := (xpart llCorner,ypart urCorner);
+ lrCorner := (xpart urCorner,ypart llCorner);
+ PageCenter := (urCorner+ulCorner+llCorner+lrCorner)/4;
+ path thebigframe;
+ thebigframe = llCorner--lrCorner--urCorner--ulCorner--cycle;
+ boolean EPSmode, Athreemode;
+ EPSmode = false;
+ Athreemode = false;
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+ def centerof( text t) =
+ begingroup
+ save momentum, counter;
+ numeric counter;
+ pair momentum;
+ counter = 0;
+ momentum = (0,0);
+ forsuffixes $=t:
+ counter := incr(counter);
+ momentum := momentum + z$;
+ endfor;
+ ( momentum/counter )
+ endgroup
enddef;
-
-% EOF
+
+ def truefyepsmode =
+ begingroup
+ EPSmode := true;
+ arestaref := arestaref/2;
+ endgroup
+ enddef;
+
+ def makearestarefathree =
+ begingroup
+ if not EPSmode:
+ Athreemode := true;
+ arestaref := arestaref*1.414;
+ fi;
+ endgroup
+ enddef;
+
+ def terminar( expr ScaleFactor, RotAngle ) =
+ begingroup
+ picture inicial,final;
+ inicial = currentpicture;
+ currentpicture := nullpicture;
+ final = inicial scaled ScaleFactor;
+ if EPSmode:
+ final := final scaled 2
+ fi;
+ final := final rotated RotAngle;
+ if Athreemode:
+ draw final shifted (PageCenter*1.414);
+ arestaref := defaultaresta;
+ Athreemode := false;
+ else:
+ draw final shifted PageCenter;
+ if EPSmode:
+ pickup grossa;
+ draw thebigframe;
+ inicial := currentpicture;
+ currentpicture := nullpicture;
+ final := inicial rotated -90;
+ draw final
+ fi;
+ fi;
+ endgroup
+ enddef;
+
+ def desvia( text zes ) =
+ begingroup
+ picture inicial,final;
+ inicial = currentpicture;
+ currentpicture := nullpicture;
+ final = inicial shifted -centerof( zes );
+ draw final
+ endgroup
+ enddef;
+
+ def desvec( expr vect ) =
+ begingroup
+ picture inicial,final;
+ inicial = currentpicture;
+ currentpicture := nullpicture;
+ final = inicial shifted -vect;
+ draw final
+ endgroup
+ enddef;
+
+% EOF \ No newline at end of file
generated by cgit v1.2.3 (git 2.39.1) at 2024-11-04 05:19:05 +0000