(FILECREATED "17-Feb-87 13:29:00" {QV}<PEDERSEN>LISP>KOTO>VECTORPLOT.;3 13149
changes to: (VARS VECTORPLOTCOMS)
previous date: "17-Feb-87 09:08:53" {QV}<PEDERSEN>LISP>KOTO>VECTORPLOT.;2)
(* Copyright (c) 1987 by Xerox Corporation. All rights reserved.)
(PRETTYCOMPRINT VECTORPLOTCOMS)
(RPAQQ VECTORPLOTCOMS ((FNS CREATEIMAGE DISPLAYIMAGE GETNEWROTATIONMATRIX MAKEDATA MAPWORLDTOWINDOW
MULTIPLYANDPROJECT THREEDPLOT THREEDPLOT-ROTATE)
(MACROS SUBMATRIXROW \MATMULT331 \MATMULT333)
(FILES CMLARRAY CMLFLOATARRAY BLAS UNBOXEDOPS)))
(DEFINEQ
(CREATEIMAGE
[LAMBDA (OLDPOINTS NEWPOINTS WINDOW N) (* jop: " 4-Feb-86 16:10")
(* *)
(TOTOPW WINDOW)
(bind (OLDBASE ←(ARRAYBASE OLDPOINTS))
(NEWBASE ←(ARRAYBASE NEWPOINTS))
(SCREEN ←(SCREENBITMAP)) for I from 1 to N as XINDEX from 0 by 2
as YINDEX from (MUL2 N) by 2
do (BITMAPBIT SCREEN (\GETBASEPTR OLDBASE XINDEX)
(\GETBASEPTR OLDBASE YINDEX)
0)
(BITMAPBIT SCREEN (\GETBASEPTR NEWBASE XINDEX)
(\GETBASEPTR NEWBASE YINDEX)
1])
(DISPLAYIMAGE
[LAMBDA (BITMAP WINDOW) (* jop: "14-Jan-85 14:22")
(* * BITMAP is an image to be displayed in WINDOW. BITMAP should be of the same dimensions as WINDOW)
(BITBLT BITMAP NIL NIL WINDOW])
(GETNEWROTATIONMATRIX
[LAMBDA (K OLDMATRIX PRODUCT) (* jop: " 3-Feb-86 20:13")
(* * Multiply OLDMATRIX by a rotation matrix K (a fixed rotation about the y axis) and return the result in
PRODUCT. K, OLDMATRIX and PRODUCT must be 3 x 3 CMLARRAYS)
(\MATMULT333 (ARRAYBASE K)
(ARRAYBASE OLDMATRIX)
(ARRAYBASE PRODUCT])
(MAKEDATA
[LAMBDA (N LOWER UPPER) (* jop: " 4-Feb-86 15:03")
(* *)
(RANDSET T)
(SETQ LOWER (FLOAT (OR LOWER 0.0)))
(SETQ UPPER (FLOAT (OR UPPER 1.0)))
(PROG [(MATRIX (MAKE-ARRAY (LIST 3 N)
(QUOTE :ELEMENT-TYPE)
(QUOTE SINGLE-FLOAT]
(for I from 0 to (SUB1 N)
do (LASET (RAND LOWER UPPER)
MATRIX 0 I)
(LASET (RAND LOWER UPPER)
MATRIX 1 I)
(LASET (RAND LOWER UPPER)
MATRIX 2 I))
(RETURN MATRIX])
(MAPWORLDTOWINDOW
[LAMBDA (PROJECTEDDATA OLDPOINTS WINDOW N) (* jop: " 4-Feb-86 17:40")
(* * IT is assumed that PROJECTEDDATA is a 2 x N ARRAY.)
(* Draw the points directly onto the screen)
(bind (SCREEN ←(SCREENBITMAP))
(PROJECTEDDATABASE ←(ARRAYBASE PROJECTEDDATA))
(OLDPBASE ←(ARRAYBASE OLDPOINTS))
X Y for I from 1 to N as XINDEX from 0 by 2 as YINDEX from (MUL2 N)
by 2
do (SETQ X (UFIX (\GETBASEFLOATP PROJECTEDDATABASE XINDEX)))
(SETQ Y (UFIX (\GETBASEFLOATP PROJECTEDDATABASE YINDEX)))
(* First erase the old point)
(BITMAPBIT SCREEN (\GETBASEPTR OLDPBASE XINDEX)
(\GETBASEPTR OLDPBASE YINDEX)
0) (* Draw the new point)
(BITMAPBIT SCREEN X Y 1) (* Save the new point)
(\PUTBASEPTR OLDPBASE XINDEX X)
(\PUTBASEPTR OLDPBASE YINDEX Y])
(MULTIPLYANDPROJECT
[LAMBDA (RMATRIX DATA AX AY PROJECTION N) (* jop: " 4-Feb-86 17:39")
(* * RMATRIX is a rotation matrix. Multiply RMATRIX and DATA then PROJECT out the Z coordinate)
(if (IGREATERP N FFTSSIZE)
then (HELP "N too large" N))
(PROG ((NLESS1 (SUB1 N))
(RMATRIXBASE (ARRAYBASE RMATRIX))
(DATABASE (ARRAYBASE DATA))
(PROJECTIONBASE (ARRAYBASE PROJECTION))
(SCRATCHBASE \SLICE1))
(* * Do the Matrix multiply)
(bind (TEMPPROJECTIONBASE ← PROJECTIONBASE)
(STARTRMATRIXOFFSET ← 0)
(TWICEN ←(MUL2 N)) for I from 0 to 1
do (bind (TEMPDATABASE ← DATABASE) for NEXTI from 0 to 2 as RMATRIXOFFSET
from STARTRMATRIXOFFSET by 2
do (* Fill SCRATCHBASE with an element from RMATRIXBASE)
(* Setup the last entry)
(\PUTBASEFLOATP SCRATCHBASE (MUL2 NLESS1)
(\GETBASEFLOATP RMATRIXBASE RMATRIXOFFSET))
(* \BLT operates backwards)
(\BLT SCRATCHBASE (\ADDBASE SCRATCHBASE 2)
(MUL2 NLESS1))
(if (EQ NEXTI 0)
then (\BLKFTIMES TEMPDATABASE SCRATCHBASE TEMPPROJECTIONBASE N)
else (\BLKFTIMES TEMPDATABASE SCRATCHBASE SCRATCHBASE N)
(\BLKFPLUS SCRATCHBASE TEMPPROJECTIONBASE TEMPPROJECTIONBASE N)
) (* increment Data row)
(SETQ TEMPDATABASE (\ADDBASE TEMPDATABASE TWICEN)))
(* increment Projecteddata row)
(SETQ TEMPPROJECTIONBASE (\ADDBASE TEMPPROJECTIONBASE TWICEN))
(* increment rotation matrix row)
(SETQ STARTRMATRIXOFFSET (\ADDBASE STARTRMATRIXOFFSET 6)))
(* * Translate)
(* fill SCRATCHBASE with AX)
(\PUTBASEFLOATP SCRATCHBASE (MUL2 NLESS1)
AX)
(\BLT SCRATCHBASE (\ADDBASE SCRATCHBASE 2)
(MUL2 NLESS1))
(\BLKFPLUS SCRATCHBASE PROJECTIONBASE PROJECTIONBASE N)
(* fill SCRATCHBASE with AY)
(\PUTBASEFLOATP SCRATCHBASE (MUL2 NLESS1)
AY)
(\BLT SCRATCHBASE (\ADDBASE SCRATCHBASE 2)
(MUL2 NLESS1)) (* now add to the second row)
(SETQ PROJECTIONBASE (\ADDBASE PROJECTIONBASE (MUL2 N)))
(\BLKFPLUS SCRATCHBASE PROJECTIONBASE PROJECTIONBASE N])
(THREEDPLOT
[LAMBDA (DATA) (* edited: "17-Feb-87 09:08")
(* * Expects DATA to be a 3 x n CMLARRAY)
(if [NOT (AND (type? ARRAY DATA)
(EQ (ARRAY-ELEMENT-TYPE DATA)
(QUOTE SINGLE-FLOAT))
(IEQP 2 (ARRAY-RANK DATA))
(IEQP 3 (ARRAY-DIMENSION DATA 0]
then (HELP "Data not in correct form"))
(* * Translate data so origin is at center of mass, and compute the maximum norm)
(PROG ((N (ARRAY-DIMENSION DATA 1))
(MAXNORM 0.0)
WINDOW SCALEDDATA MX AX MY AY)
(SETQ SCALEDDATA (MAKE-ARRAY (LIST 3 N)
(QUOTE :ELEMENT-TYPE)
(QUOTE SINGLE-FLOAT)))
(* * Normalize to zero mean)
(* Copy DATA to SCALEDDATA)
(BLAS.ARRAYBLT DATA 0 1 SCALEDDATA 0 1)
(BLAS.ADD (FMINUS (FQUOTIENT (BLAS.SUM SCALEDDATA 0 1 N)
N))
SCALEDDATA 0 1 N)
(BLAS.ADD (FMINUS (FQUOTIENT (BLAS.SUM SCALEDDATA N 1 N)
N))
SCALEDDATA N 1 N)
(BLAS.ADD (FMINUS (FQUOTIENT (BLAS.SUM SCALEDDATA (ITIMES 2 N)
1 N)
N))
SCALEDDATA
(ITIMES 2 N)
1 N) (* find the largest two norm)
[SETQ MAXNORM (bind (MAX ← MIN.FLOAT)
DOTP for I from 0 to (SUB1 N)
do (SETQ DOTP (BLAS.DOTPROD SCALEDDATA I N SCALEDDATA I N 3))
(if (FGREATERP DOTP MAX)
then (SETQ MAX DOTP))
finally (RETURN (SQRT MAX]
(* * Compute a WORLD to WINDOW transformation so that twice MAXNORM will always fit in the window)
(SETQ WINDOW (CREATEW NIL "Three D plot"))
[LET ((WWIDTH (WINDOWPROP WINDOW (QUOTE WIDTH)))
(WHEIGHT (WINDOWPROP WINDOW (QUOTE HEIGHT)))
SIZE SCALE TRANSLATE) (* make data fit in largest square sub region)
(SETQ SIZE (IMIN WWIDTH WHEIGHT))
(SETQ SCALE (QUOTIENT SIZE (TIMES 2.0 MAXNORM)))
(* Prescale the data)
(BLAS.SCAL SCALE SCALEDDATA 0 1 N)
(BLAS.SCAL SCALE SCALEDDATA N 1 N)
(BLAS.SCAL SCALE SCALEDDATA (ITIMES 2 N)
1 N)
(SETQ TRANSLATE (QUOTIENT SIZE 2.0))
(SETQ AX (PLUS TRANSLATE (QUOTIENT (DIFFERENCE WWIDTH SIZE)
2.0)))
(SETQ AY (PLUS TRANSLATE (QUOTIENT (DIFFERENCE WHEIGHT SIZE)
2.0]
(THREEDPLOT-ROTATE SCALEDDATA AX AY WINDOW])
(THREEDPLOT-ROTATE
[LAMBDA (DATA AX AY WINDOW) (* jop: " 4-Feb-86 17:40")
(* *)
(PROG ((N (ARRAY-DIMENSION DATA 1))
[OLDROTATIONMATRIX (MAKE-ARRAY (QUOTE (3 3))
(QUOTE :ELEMENT-TYPE)
(QUOTE SINGLE-FLOAT)
(QUOTE :INITIAL-CONTENTS)
(QUOTE ((1.0 0.0 0.0)
(0.0 1.0 0.0)
(0.0 0.0 1.0]
[NEWROTATIONMATRIX (MAKE-ARRAY (QUOTE (3 3))
(QUOTE :ELEMENT-TYPE)
(QUOTE SINGLE-FLOAT)
(QUOTE :INITIAL-CONTENTS)
(QUOTE ((1.0 0.0 0.0)
(0.0 1.0 0.0)
(0.0 0.0 1.0]
[CONSTANTYAXIS (MAKE-ARRAY (QUOTE (3 3))
(QUOTE :ELEMENT-TYPE)
(QUOTE SINGLE-FLOAT)
(QUOTE :INITIAL-CONTENTS)
(LIST (LIST (COS 2.5)
0.0
(FMINUS (SIN 2.5)))
(LIST 0.0 1.0 0.0)
(LIST (SIN 2.5)
0.0
(COS 2.5]
[CONSTANTXAXIS (MAKE-ARRAY (QUOTE (3 3))
(QUOTE :ELEMENT-TYPE)
(QUOTE SINGLE-FLOAT)
(QUOTE :INITIAL-CONTENTS)
(LIST (LIST 1.0 0.0 0.0)
(LIST 0.0 (COS 2.5)
(SIN 2.5))
(LIST 0.0 (FMINUS (SIN 2.5))
(COS 2.5]
PROJECTEDDATA OLDPOINTS)
(SETQ PROJECTEDDATA (MAKE-ARRAY (LIST 2 N)
(QUOTE :ELEMENT-TYPE)
(QUOTE SINGLE-FLOAT)))
(SETQ OLDPOINTS (MAKE-ARRAY (LIST 2 N)))
(* * Adjust for position of WINDOW)
(SETQ AX (PLUS AX (DSPXOFFSET NIL WINDOW)))
(SETQ AY (PLUS AY (DSPYOFFSET NIL WINDOW)))
(* * Initialize OLDPOINTS)
(bind (IAX ←(FIXR AX))
(IAY ←(FIXR AY)) for I from 0 to (SUB1 N)
do (PASET IAX OLDPOINTS 0 I)
(PASET IAY OLDPOINTS 1 I))
(* * Initialize the display)
(MULTIPLYANDPROJECT OLDROTATIONMATRIX DATA AX AY PROJECTEDDATA N)
(MAPWORLDTOWINDOW PROJECTEDDATA OLDPOINTS WINDOW N)
(* * go into an loop, refreshing the display at each iteration)
(bind (OLDMATRIX ← OLDROTATIONMATRIX)
(NEWMATRIX ← NEWROTATIONMATRIX)
(K ← CONSTANTYAXIS)
(BEFORE ←(CLOCK 0))
(NUMITER ← 0)
FLIP TTIME ENDITER while (NOT ENDITER)
do (if (MOUSESTATE LEFT)
then (SETQ K CONSTANTYAXIS)
elseif (MOUSESTATE RIGHT)
then (SETQ K CONSTANTXAXIS)
elseif (MOUSESTATE MIDDLE)
then (SETQ ENDITER T))
(SETQ NUMITER (ADD1 NUMITER))
(GETNEWROTATIONMATRIX K OLDMATRIX NEWMATRIX)
(MULTIPLYANDPROJECT NEWMATRIX DATA AX AY PROJECTEDDATA N)
(MAPWORLDTOWINDOW PROJECTEDDATA OLDPOINTS WINDOW N)
(if (SETQ FLIP (NOT FLIP))
then (SETQ OLDMATRIX NEWROTATIONMATRIX)
(SETQ NEWMATRIX OLDROTATIONMATRIX)
else (SETQ OLDMATRIX OLDROTATIONMATRIX)
(SETQ NEWMATRIX NEWROTATIONMATRIX))
finally (SETQ TTIME (FQUOTIENT (DIFFERENCE (CLOCK 0)
BEFORE)
1000))
(PRINTOUT T N " Points at " (FQUOTIENT NUMITER TTIME)
, "Iterations per second" T])
)
(DECLARE: EVAL@COMPILE
(PUTPROPS SUBMATRIXROW MACRO (OPENLAMBDA (MATRIX MATRIXROW SCALAR RESULT RESULTROW NELTS)
(ADDMATRIXROW MATRIX MATRIXROW (FMINUS SCALAR)
RESULT RESULTROW NELTS)))
(PUTPROPS \MATMULT331 DMACRO ((A B C)
((OPCODES UBFLOAT3 4)
A B C)))
(PUTPROPS \MATMULT333 DMACRO ((A B C)
((OPCODES UBFLOAT3 1)
A B C)))
)
(FILESLOAD CMLARRAY CMLFLOATARRAY BLAS UNBOXEDOPS)
(PUTPROPS VECTORPLOT COPYRIGHT ("Xerox Corporation" 1987))
(DECLARE: DONTCOPY
(FILEMAP (NIL (578 12653 (CREATEIMAGE 588 . 1200) (DISPLAYIMAGE 1202 . 1470) (GETNEWROTATIONMATRIX
1472 . 1873) (MAKEDATA 1875 . 2491) (MAPWORLDTOWINDOW 2493 . 3639) (MULTIPLYANDPROJECT 3641 . 6460) (
THREEDPLOT 6462 . 9269) (THREEDPLOT-ROTATE 9271 . 12651)))))
STOP