summaryrefslogtreecommitdiff
path: root/Master/texmf-dist/dvips/pst-3dplot/pst-3dplot.pro
diff options
context:
space:
mode:
authorKarl Berry <karl@freefriends.org>2007-12-23 01:05:36 +0000
committerKarl Berry <karl@freefriends.org>2007-12-23 01:05:36 +0000
commit2154eaa8df555b17f72558d748c04bfe1c9c144b (patch)
tree6db55c6713bdc0b1b6428a94eeaa969c45042c85 /Master/texmf-dist/dvips/pst-3dplot/pst-3dplot.pro
parentf0572c6c196ee131830234fbde3a305ecdced4f0 (diff)
pst-3dplot 1.78 (23dec07)
git-svn-id: svn://tug.org/texlive/trunk@5827 c570f23f-e606-0410-a88d-b1316a301751
Diffstat (limited to 'Master/texmf-dist/dvips/pst-3dplot/pst-3dplot.pro')
-rw-r--r--Master/texmf-dist/dvips/pst-3dplot/pst-3dplot.pro381
1 files changed, 321 insertions, 60 deletions
diff --git a/Master/texmf-dist/dvips/pst-3dplot/pst-3dplot.pro b/Master/texmf-dist/dvips/pst-3dplot/pst-3dplot.pro
index 12f52ad1443..251f0612ff5 100644
--- a/Master/texmf-dist/dvips/pst-3dplot/pst-3dplot.pro
+++ b/Master/texmf-dist/dvips/pst-3dplot/pst-3dplot.pro
@@ -15,22 +15,24 @@
%% `pst-3dplot' is a PSTricks package to draw 3d curves and graphical objects
%%
%%
-%% version 0.23 / 2007-08-18 Herbert Voss <voss _at_ PSTricks.de>
+%% version 0.25 / 2007-12-22 Herbert Voss <voss _at_ PSTricks.de>
+%% with contributions of Darrell Lamm <darrell.lamm _at_ gtri.gatech.edu<
+%%
%
-/tx@3DPlotDict 50 dict def
+/tx@3DPlotDict 200 dict def
tx@3DPlotDict begin
%
/printDot { gsave 2 copy 2 0 360 arc fill stroke grestore } def
%
/saveCoor {
- dzUnit mul /z exch def
- dyUnit mul /y exch def
- dxUnit mul /x exch def
+ dzUnit mul /z ED
+ dyUnit mul /y ED
+ dxUnit mul /x ED
} def
%
/ConvertTo2D {
RotatePoint
- 1 { % dummy loop, will run only 1 time
+ 1 { % dummy loop, will run only 1 time, allows exit
coorType 0 le { % the default |
/x2D x leftHanded not { neg } if Alpha cos mul y Alpha sin mul add def % /\ co system
/y2D x leftHanded { neg } if Alpha sin mul y Alpha cos mul add neg Beta sin mul z Beta cos mul add def
@@ -43,6 +45,10 @@ tx@3DPlotDict begin
/x2D y x 0.5 mul sub def
/y2D z x 0.5 mul sub def
exit } if
+ coorType 3 le { % coorType |/_ with a 1/sqrt(2) shortend x-axis and 135 degrees
+ /x2D y x -0.5 mul sub def
+ /y2D z x -0.5 mul sub def
+ exit } if
} repeat
} def
%
@@ -75,7 +81,7 @@ tx@3DPlotDict begin
%
/convertStackTo2D {
counttomark
- /n exch def /n3 n 3 div cvi def
+ /n ED /n3 n 3 div cvi def
n3 {
n -3 roll
SphericalCoor { ConvertToCartesian } { saveCoor } ifelse
@@ -85,43 +91,181 @@ tx@3DPlotDict begin
} repeat
} def
%
+/RotSet (set ) def
+%
+/eulerRotation false def
+% Matrix multiplication procedure
+/matmul {
+
+ /M@tMulDict 20 dict def
+ M@tMulDict begin
+ /m2 ED
+ /m1 ED
+ m1 dup length 2 sub 2 getinterval aload pop
+ /col1max ED
+ /row1max ED
+ m2 dup length 2 sub 2 getinterval aload pop
+ /col2max ED
+ /row2max ED
+ /m3 row1max col2max mul 2 add array def
+ m3 dup length 2 sub row1max col2max 2 array astore putinterval
+ 0 1 row1max 1 sub {
+ /row ED
+ 0 1 col2max 1 sub {
+ /col ED
+ /sum 0 def
+ 0 1 col1max 1 sub{
+ /rowcol ED
+ sum
+ m1 row col1max mul rowcol add get
+ m2 rowcol col2max mul col add get
+ mul add
+ /sum ED
+ } for
+ m3 row col2max mul col add sum put
+ } for
+ } for
+ m3
+ end % end of M@tMulDict
+
+} def
+%
+/SetMQuaternion {
+
+ /MnewTOold 11 array def
+
+ /Qu@ternionDict 30 dict def
+ Qu@ternionDict begin
+
+ /normRotVec xRotVec yRotVec zRotVec 3 array astore VecNorm def
+ normRotVec 0 gt
+ {/xRotVecNorm xRotVec normRotVec div def
+ /yRotVecNorm yRotVec normRotVec div def
+ /zRotVecNorm zRotVec normRotVec div def
+ RotAngle}
+ {/xRotVecNorm 1 def
+ /yRotVecNorm 0 def
+ /zRotVecNorm 0 def
+ 0} ifelse
+
+ 2 div dup
+ /q0 exch cos def
+ sin dup dup
+ /q1 exch xRotVecNorm mul def
+ /q2 exch yRotVecNorm mul def
+ /q3 exch zRotVecNorm mul def
+
+ /q0q0 q0 q0 mul def
+ /q0q1 q0 q1 mul def
+ /q0q2 q0 q2 mul def
+ /q0q3 q0 q3 mul def
+
+ /q1q1 q1 q1 mul def
+ /q1q2 q1 q2 mul def
+ /q1q3 q1 q3 mul def
+
+ /q2q2 q2 q2 mul def
+ /q2q3 q2 q3 mul def
+
+ /q3q3 q3 q3 mul def
+
+ MnewTOold 0 q0q0 q1q1 add q2q2 sub q3q3 sub put
+ MnewTOold 1 q1q2 q0q3 sub 2 mul put
+ MnewTOold 2 q1q3 q0q2 add 2 mul put
+
+ MnewTOold 3 q1q2 q0q3 add 2 mul put
+ MnewTOold 4 q0q0 q1q1 sub q2q2 add q3q3 sub put
+ MnewTOold 5 q2q3 q0q1 sub 2 mul put
+
+ MnewTOold 6 q1q3 q0q2 sub 2 mul put
+ MnewTOold 7 q2q3 q0q1 add 2 mul put
+ MnewTOold 8 q0q0 q1q1 sub q2q2 sub q3q3 add put
+
+ MnewTOold 9 3 put
+ MnewTOold 10 3 put
+
+ end % end of Qu@ternionDict
+
+} def
+%
+/SetMxyz {
+
+ 1.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 1.0 3 3 11 array astore /MnewTOold ED
+ RotSequence cvx exec % Now create a new MnewTOold using xyz, etc.
+
+} def
%
-/CalcCoordinates {% from pst-vue3d
- formulesTroisD
- Xi 28.45 mul Yi 28.45 mul
+/ConcatMQuaternion {
+
+ MnewTOold % Push onto stack
+ SetMQuaternion % Uses [xyz]RotVec and RotAngle to make MnewToOld
+ MnewTOold matmul /MnewTOold ED
+
+} def
+%
+/ConcatMxyz {
+
+ MnewTOold % Push onto stack
+ SetMxyz % Uses RotX, etc. to set MnewTOold
+ MnewTOold matmul /MnewTOold ED
+
} def
-% pour la 3D conventionnelle
-/formulesTroisD{%
- /xObservateur x Sin1 mul neg y Cos1 mul add def
- /yObservateur x Cos1Sin2 mul neg y Sin1Sin2 mul sub z Cos2 mul add def
- /zObservateur x neg Cos1Cos2 mul y Sin1Cos2 mul sub z Sin2 mul sub def
- /Xi xObservateur zObservateur Div def
- /Yi yObservateur zObservateur Div def
+%
+/RotatePoint{
+ MnewTOold x y z 3 1 5 array astore matmul
+ 0 3 getinterval aload pop
+ /z ED
+ /y ED
+ /x ED
} def
%
-/RotatePointXYZ{% Mxx are defined in the TeX file
- /xi M11 x mul M12 y mul add M13 z mul add def
- /yi M21 x mul M22 y mul add M23 z mul add def
- /zi M31 x mul M32 y mul add M33 z mul add def
- /x xi def
- /y yi def
- /z zi def
+/makeMoldTOnew {
+ /MoldTOnew 11 array def
+ MoldTOnew 0 MnewTOold 0 get put
+ MoldTOnew 1 MnewTOold 3 get put
+ MoldTOnew 2 MnewTOold 6 get put
+ MoldTOnew 3 MnewTOold 1 get put
+ MoldTOnew 4 MnewTOold 4 get put
+ MoldTOnew 5 MnewTOold 7 get put
+ MoldTOnew 6 MnewTOold 2 get put
+ MoldTOnew 7 MnewTOold 5 get put
+ MoldTOnew 8 MnewTOold 8 get put
+ MoldTOnew 9 3 put
+ MoldTOnew 10 3 put
} def
%
/RotXaxis {
- /yTemp y RotX cos mul z RotX sin mul sub def
- /z y RotX sin mul z RotX cos mul add def
- /y yTemp def
+ eulerRotation
+ {1 0 0}
+ {makeMoldTOnew MoldTOnew 1 0 0 3 1 5 array astore matmul
+ 0 3 getinterval aload pop} ifelse
+ /zRotVec ED
+ /yRotVec ED
+ /xRotVec ED
+ /RotAngle RotX def
+ ConcatMQuaternion
} def
/RotYaxis {
- /xTemp x RotY cos mul z RotY sin mul add def
- /z x RotY sin mul neg z RotY cos mul add def
- /x xTemp def
+ eulerRotation
+ {0 1 0}
+ {makeMoldTOnew MoldTOnew 0 1 0 3 1 5 array astore matmul
+ 0 3 getinterval aload pop} ifelse
+ /zRotVec ED
+ /yRotVec ED
+ /xRotVec ED
+ /RotAngle RotY def
+ ConcatMQuaternion
} def
/RotZaxis {
- /xTemp x RotZ cos mul y RotZ sin mul sub def
- /y x RotZ sin mul y RotZ cos mul add def
- /x xTemp def
+ eulerRotation
+ {0 0 1}
+ {makeMoldTOnew MoldTOnew 0 0 1 3 1 5 array astore matmul
+ 0 3 getinterval aload pop} ifelse
+ /zRotVec ED
+ /yRotVec ED
+ /xRotVec ED
+ /RotAngle RotZ def
+ ConcatMQuaternion
} def
/xyz { RotXaxis RotYaxis RotZaxis } def
/yxz { RotYaxis RotXaxis RotZaxis } def
@@ -129,13 +273,13 @@ tx@3DPlotDict begin
/xzy { RotXaxis RotZaxis RotYaxis } def
/zxy { RotZaxis RotXaxis RotYaxis } def
/zyx { RotZaxis RotYaxis RotXaxis } def
-%
-/RotatePoint { RotSequence cvx exec } def
+/quaternion { } def % Null
%
/VecNorm { 0 exch { dup mul add } forall sqrt } def
%
/UnitVec { % on stack is [a]; returns a vector with [a][a]/|a|=1
dup VecNorm /norm ED
+ norm 0 lt {/norm 0 def} if
{ norm div } forall 3 array astore } def
%
/AxB { % on the stack are the two vectors [a][b]
@@ -169,10 +313,10 @@ tx@3DPlotDict begin
% les rayons de lumière
xLight dup mul yLight dup mul zLight dup mul add add sqrt /NormeLight ED
% the color values
- /K exch def
- /Yellow exch def
- /Magenta exch def
- /Cyan exch def
+ /K ED
+ /Yellow ED
+ /Magenta ED
+ /Cyan ED
} def
%
/facetteSphere {
@@ -182,28 +326,28 @@ tx@3DPlotDict begin
/Zpoint Rsphere phi sin mul CZ add def
Xpoint Ypoint Zpoint tx@ProjThreeD moveto
theta 1 theta increment add {%
- /theta1 exch def
+ /theta1 ED
/Xpoint Rsphere theta1 cos mul phi cos mul CX add def
/Ypoint Rsphere theta1 sin mul phi cos mul CY add def
/Zpoint Rsphere phi sin mul CZ add def
Xpoint Ypoint Zpoint tx@ProjThreeD lineto
} for
phi 1 phi increment add {
- /phi1 exch def
+ /phi1 ED
/Xpoint Rsphere theta increment add cos mul phi1 cos mul CX add def
/Ypoint Rsphere theta increment add sin mul phi1 cos mul CY add def
/Zpoint Rsphere phi1 sin mul CZ add def
Xpoint Ypoint Zpoint tx@ProjThreeD lineto
} for
theta increment add -1 theta {%
- /theta1 exch def
+ /theta1 ED
/Xpoint Rsphere theta1 cos mul phi increment add cos mul CX add def
/Ypoint Rsphere theta1 sin mul phi increment add cos mul CY add def
/Zpoint Rsphere phi increment add sin mul CZ add def
Xpoint Ypoint Zpoint tx@ProjThreeD lineto
} for
phi increment add -1 phi {
- /phi1 exch def
+ /phi1 ED
/Xpoint Rsphere theta cos mul phi1 cos mul CX add def
/Ypoint Rsphere theta sin mul phi1 cos mul CY add def
/Zpoint Rsphere phi1 sin mul CZ add def
@@ -215,17 +359,17 @@ tx@3DPlotDict begin
/MaillageSphere {
% on stack must be x y z Radius increment C M Y K
setColorLight
- /increment exch def
- /Rsphere exch def
- /CZ exch def
- /CY exch def
- /CX exch def
+ /increment ED
+ /Rsphere ED
+ /CZ ED
+ /CY ED
+ /CX ED
/StartTheta 0 def
/condition { PSfacetteSphere 0 ge } def
-90 increment 90 increment sub {%
- /phi exch def
+ /phi ED
StartTheta increment 360 StartTheta add increment sub {%
- /theta exch def
+ /theta ED
% Centre de la facette
/Xpoint Rsphere theta increment 2 div add cos mul phi increment 2 div add cos mul CX add def
/Ypoint Rsphere theta increment 2 div add sin mul phi increment 2 div add cos mul CY add def
@@ -267,7 +411,7 @@ tx@3DPlotDict begin
/PlanCoupeCylinder { %
/TableauxPoints [
0 1 359 {
- /phi exch def
+ /phi ED
[ Radius phi Height ConvCyl2d ] % on décrit le cercle
} for
] def
@@ -374,7 +518,7 @@ tx@3DPlotDict begin
/height ED
/radius ED
startAngle increment endAngle {
- /Angle exch def
+ /Angle ED
radius Angle 0 ConvCylToCartesian MoveTo
radius Angle height ConvCylToCartesian LineTo
} for
@@ -494,13 +638,13 @@ tx@3DPlotDict begin
/V {Z sqrt} bind def
/TableauxPoints [
0 1 359 {
- /U exch def [ U U Z V calculate2DPoint ] % on décrit le cercle
+ /U ED [ U U Z V calculate2DPoint ] % on décrit le cercle
} for
] def
newpath
TableauxPoints 0 get aload pop moveto
0 1 359 {
- /compteur exch def
+ /compteur ED
TableauxPoints compteur get aload pop
lineto } for
closepath
@@ -510,20 +654,20 @@ tx@3DPlotDict begin
newpath
U U Z V calculate2DPoint moveto
U 1 U increment add {%
- /U1 exch def
+ /U1 ED
U1 U1 Z V calculate2DPoint lineto
} for
Z pas10 Z pas add pas10 add{
- /Z1 exch def
+ /Z1 ED
/V {Z1 sqrt} bind def
U1 U1 Z1 V calculate2DPoint lineto
} for
U increment add -1 U {%
- /U2 exch def
+ /U2 ED
U2 U2 Z pas add V calculate2DPoint lineto
} for
Z pas add pas10 sub pas10 neg Z pas10 sub {
- /Z2 exch def
+ /Z2 ED
/V Z2 abs sqrt def
U U Z2 V calculate2DPoint lineto
} for
@@ -539,7 +683,7 @@ tx@3DPlotDict begin
/Z ED
/V Z sqrt def
0 increment 360 increment sub {%
- /U exch def
+ /U ED
% Centre de la facette
/Ucentre U increment 2 div add def
/Vcentre Z pas 2 div add sqrt def
@@ -576,5 +720,122 @@ tx@3DPlotDict begin
} for
} def
%
+% ------------------------------------ math stuff ----------------------------------
+%
+% Matrix A in arrays of rows A[[row1][row2]...]
+% with [row1]=[a11 a12 ... b1]
+% returns on stack solution vector X=[x1 x2 ... xn]
+/SolveLinEqSystem { % on stack matrix M=[A,b] (A*x=b)
+ 10 dict begin % hold all ocal
+ /A exch def
+ /Rows A length def % Rows = number of rows
+ /Cols A 0 get length def % Cols = number of columns
+ /Index [ 0 1 Rows 1 sub { } for ] def % Index = [0 1 2 ... Rows-1]
+ /col 0 def
+ /row 0 def
+ /PR Rows array def % PR[c] = pivot row for row row
+ { % starts the loop, find pivot entry in row r
+ col Cols ge row Rows ge or { exit } if % col < Cols and row < Rows else exit
+ /pRow row def % pRow = pivot row
+ /max A row get col get abs def % get A[row[col]], first A[0,0]
+ row 1 add 1 Rows 1 sub { % starts for loop 1 1 Rows-1
+ /j exch def % index counter
+ /x A j get col get abs def % get A[j[r]]
+ x max gt { % x>max, then save position
+ /pRow j def
+ /max x def
+ } if
+ } for % now we have the row with biggest A[0,1]
+ % with pRow = the pivot row
+ max 0 gt { % swap entries pRow and row in i
+ /tmp Index row get def
+ Index row Index pRow get put
+ Index pRow tmp put % and columns pRow and row in A
+ /tmp A row get def
+ A row A pRow get put
+ A pRow tmp put % pivot
+ /row0 A row get def % the pivoting row
+ /p0 row0 col get def % the pivot value
+ row 1 add 1 Rows 1 sub { % start for loop
+ /j exch def
+ /c1 A j get def
+ /p c1 col get p0 div def
+ c1 col p put % subtract (p1/p0)*row[i] from row[j]
+ col 1 add 1 Cols 1 sub { % start for loop
+ /i exch def
+ c1 dup i exch % c1 i c1
+ i get row0 i get p mul sub put
+ } for
+ } for
+ PR row col put
+ /col col 1 add def
+ /row row 1 add def
+ }{ % all zero entries
+ /row row 1 add def % continue loop with same row
+ } ifelse
+ } loop
+ /X A def % solution vector
+ A Rows 1 sub get dup
+ Cols 1 sub get exch
+ Cols 2 sub get div
+ X Rows 1 sub 3 -1 roll put % X[n]
+ Rows 2 sub -1 0 { % for loop to calculate X[i]
+ /xi exch def % current index
+ A xi get % i-th row
+ /Axi exch def
+ /sum 0 def
+ Cols 2 sub -1 xi 1 add {
+ /n exch def
+ /sum sum Axi n get X n get mul add def
+ } for
+ Axi Cols 1 sub get % b=Axi[Cols-1]
+ sum sub % b-sum
+ Axi xi get div % b-sum / Axi[xi]
+ X xi 3 -1 roll put % X[xi]
+ } for
+ X
+ end
+} def
%
-end
+% u v -> u+v
+/vector-add { 1 dict begin
+ /v exch def
+ [ exch
+ 0 % u i
+ exch { % i u[i]
+ v % i u[i] v
+ 2 index get add % i u[i]+v[i]
+ exch 1 add % i
+ } forall
+ pop
+ ]
+ end
+} def
+
+% u v -> u-v
+/vector-sub { 1 dict begin
+ /v exch def
+ [ exch
+ 0 % u i
+ exch { % i u[i]
+ v % i u[i] v
+ 2 index get sub % i u[i]+v[i]
+ exch 1 add % i
+ } forall
+ pop
+ ]
+end } def
+%
+% [v] c -> [c.v]
+/vector-scale { 1 dict begin
+ /c exch def
+ [ exch
+ { % s i u[i]
+ c mul % s i u[i] v
+ } forall
+ ]
+ end } def
+%
+%
+end % tx@3DPlotDict
+