C ALGORITHM 483 COLLECTED ALGORITHMS FROM ACM. C ALGORITHM APPEARED IN COMM. ACM, VOL. 17, NO. 09, C P. 520. C THIS PROGRAM GENERATES AN EXAMPLE OF A CROSSHATCHED A 10 C FIGURE, THAT IS, ONE FIGURE WHOSE LINES RUN PARALLEL TO A 20 C THE X-AXIS OVERLAYED BY ANOTHER FIGURE WHOSE LINES RUN A 30 C PARALLEL TO THE Z-AXIS. THE FUNCTION IS A PRODUCT TO TWO A 40 C SINC (I.E. SINF(X)/X) FUNCTIONS. A 50 DIMENSION MASK(2000), VERTEX(16), OUTBUF(61), Z(61) A 60 C FIRST FIGURE A 70 C GENERATE DATA RUNNING PARALLEL TO X-AXIS A 80 DO 20 NLINE=1,61 A 90 BEAMV = SINC(15.0*SINF((3*NLINE-93)*0.017453293)) A 100 DO 10 NPOINT=1,61 A 110 OUTBUF(NPOINT) = A 120 * BEAMV*SINC(7.5*SINF((3*NPOINT-93)*0.017453293)) + A 130 * 0.25 A 140 10 CONTINUE A 150 C PLOT EACH LINE AS IT IS COMPUTED A 160 CALL PLOT3D(10, 0.0, OUTBUF, 0.0, 0.1, 4.0, -0.1, A 170 * NLINE, 61, -45., -45., 5.0, 3.0, 10.0, MASK, 0) A 180 20 CONTINUE A 190 C SECOND FIGURE A 200 C GENERATE ARRAY OF Z-COMPONENTS A 210 DO 30 NLINE=1,61 A 220 Z(NLINE) = -0.1*(NLINE-1) A 230 30 CONTINUE A 240 C GENERATE DATA RUNNING PARALLEL TO Z-AXIS A 250 DO 50 NLINE=1,61 A 260 X = 0.1*(NLINE-1) A 270 BEAMV = SINC(7.5*SINF((3*NLINE-93)*0.017453293)) A 280 DO 40 NPOINT=1,61 A 290 OUTBUF(NPOINT) = A 300 * BEAMV*SINC(15.0*SINF((3*NPOINT-93)*0.017453293)) + A 310 * 0.25 A 320 40 CONTINUE A 330 C PLOT EACH LINE AS IT IS COMPUTED A 340 CALL PLOT3D(1011, X, OUTBUF, Z, 0.0, 4.0, 1.0, A 350 * NLINE, 61, -45., -45., 5.0, 3.0, 10.0, MASK, VERTEX) A 360 50 CONTINUE A 370 C DRAW A FRAME ON THE FIGURE A 380 CALL FRAMER(3, VERTEX, MASK) A 390 STOP A 400 END A 410 SUBROUTINE PLOT3D(IVXYZ, XDATA, YDATA, ZDATA, XSCALE, B 10 * YSCALE, ZSCALE, NLINE, NPNTS, PHI, THETA, XREF, * YREF, XLENTH, MASK, VERTEX) C MASKED 3-DIMENSIONAL PLOT PROGRAM WITH ROTATIONS C THIS ROUTINE WILL ACCEPT 3-DIMENSIONAL DATA IN VARIOUS C FORMS AS INPUT, ROTATE IT IN 3-SPACE TO ANY ANGLE, C AND PLOT THE PROJECTION OF THE RESULTING FIGURE ONTO THE C XY PLANE. LINEAR INTERPOLATION IS USED BETWEEN DATA C POINTS. THOSE LINES OF A FIGURE WHICH SHOULD BE HIDDEN BY C A PREVIOUS LINE ARE MASKED. C THE MASKING TECHNIQUE USED BY THIS ROUTINE IS BASED ON C TWO PREMISES - C LINES IN THE FOREGROUND (POSITIVE Z DIRECTION) C ARE PLOTTED BEFORE LINES IN THE BACKGROUND. C A LINE OR PORTION OF A LINE IS MASKED (HIDDEN) IF C IT LIES WITHIN THE REGION BOUNDED BY PREVIOUSLY C PLOTTED LINES. C EACH CALL TO PLOT3D CAUSES ONE LINE OF A FIGURE TO BE C PLOTTED. C TWO PARAMETERS OF THE PLOTTER ARE SET ON THE INITIAL CALL C FOR EACH FIGURE - C (PIPI) IS THE NUMBER OF PLOTTER INCREMENTS PER INCH. C (NYPI) IS THE NUMBER OF INCREMENTS AVAILABLE ACROSS THE C WIDTH OF THE PAPER (Y-DIRECTION). C WHEN A NEW FIGURE IS INITIATED, THE PLOTTER ORIGIN IS SET C AT THE BOTTOM OF THE PAPER BY PLOT3D AND SHOULD NOT BE C MOVED UNTIL THE FIGURE IS COMPLETED. C INPUT PARAMETERS - C (IVXYZ) IS A FOUR DIGIT DECIMAL INTEGER WHICH IS USED TO C SELECT VARIOUS INPUT/OUTPUT OPTIONS. THESE DIGITS, IN C DECREASING ORDER OF MAGNITUDE, WILL BE REFERRED TO AS V, C X, Y, AND Z. C IF V .NE. 0, THE VERTICES OF THE CURRENT FIGURE AND THEIR C PROJECTION ONTO THE Y=0 PLANE, WILL BE STORED IN A 16 C ENTRY REAL ARRAY (VERTEX), AND WILL BE UPDATED AS EACH C LINE IS PLOTTED. THESE COORDINATES ARE IN INCHES AND C RELATIVE TO THE CURRENT PLOTTER ORIGIN. THE X Y PAIRS C ARE ORDERED SO THAT THE FIRST PAIR COORESPONDS TO THE C FIRST POINT OF THE FIGURE, THE SECOND PAIR COORESPONDS C TO THE LAST POINT OF THE FIRST LINE, AND THE FOLLOWING C PAIRS ARE ORDERED IN A CIRCULAR FASHION. THE PAIRS ON THE C Y=0 PLANE OF THE FIGURE, THEN FOLLOW IN THE SAME ORDER. C IF V=0, THE VERTEX PARAMETER IS IGNORED, BUT SHOULD NOT C BE DELETED C IF X=0, THE X-COMPONENTS OF THIS LINE ARE ASSUMED TO BE C EQUALLY SPACED, AND ARE COMPUTED BY C X(I)=XDATA+(I-1)*XSCALE C WHERE (XDATA) IS THE INITIAL VALUE IN INCHES AND (XSCALE) C IS THE SPACING BETWEEN POINTS IN INCHES. IF X .NE. 0, THE C X-COMPONENTS OF THIS LINE ARE READ FROM AN ARRAY AND C MODIFIED BY C X(I)=XDATA(I)*XSCALE C WHERE (XSCALE) IS A SCALE FACTOR. C THE SAME RELATIONS HOLD FOR THE Y-COMPONENTS, THAT IS, IF C Y=0 C Y(I)=YDATA+(I-1)*YSCALE C AND IF Y .NE. 0 C Y(I)=YDATA(I)*YSCALE C IF Z=0, THE Z-COMPONENTS OF THIS LINE ARE ALL ASSUMED TO C BE EQUAL, AND ARE COMPUTED BY C Z(I)=ZDATA+(NLINE-1)*ZSCALE C WHERE (NLINE) IS SOME INTEGER ASSOCIATED WITH THIS LINE. C IF Z .NE. 0, AGAIN WE HAVE C Z(I)=ZDATA(I)*ZSCALE C WHEN (NLINE) IS EQUAL TO ONE, IT INDICATES THE BEGINNING C OF A NEW FIGURE. A CALL TO PLOT3D WITH (NLINE) EQUAL TO C ZERO BEFORE INITIATING A NEW FIGURE SIMULATES A LINE DRAWN C AT THE BOTTOM OF THE PAGE. THEREFORE ONLY THOSE PORTIONS C OF A LINE LYING ABOVE ALL PREVIOUS LINES WILL BE PLOTTED. C ALL OTHER PARAMETERS ARE IGNORED ON SUCH A CALL. C (NPNTS) IS THE NUMBER OF POINTS ON THIS LINE, AND MAY BE C ALTERED FROM LINE TO LINE. C (PHI) AND (THETA) ARE THE TWO ANGLES (IN DEGREES) USED TO C SPECIFY THE DESIRED 3-DIMENSIONAL ROTATION. THE FOLLOWING C TWO DEFINIATIONS OF THESE ROTATIONS ARE EQUIVALENT - C IN TERMS OF ROTATIONS OF AXES, THE INITIAL SYSTEM OF AXES, C XYZ, IS ROTATED BY AN ANGLE (PHI) COUNTERCLOCKWISE ABOUT C THE Y-AXIS, AND THE RESULTANT SYSTEM IS LABELED THE TUV C AXES. THE TUV AXES ARE THEN ROTATED BY AN ANGLE (THETA) C COUNTERCLOCKWISE ABOUT THE T-AXIS, AND THIS FINAL SYSTEM C IS LABELED THE PQR AXES. THE PLOTTED FIGURE IS THE C PROJECTION OF THE ORIGINAL FIGURE ONTO THE PQ-PLANE. C IN TERMS OF ROTATIONS OF COORDINATES, THE FIGURE IS FIRST C ROTATED BY AN ANGLE (THETA) CLOCKWISE ABOUT THE X-AXIS. C THE RESULTANT FIGURE IS THEN ROTATED BY AN ANGLE (PHI) C CLOCKWISE ABOUT ITS Y-AXIS. THE PLOTTED FIGURE IS THE C PROJECTION OF THIS FINAL FIGURE ONTO THE XY-PLANE. C WARNING. SOME ROTATIONS WILL ALTER THE FOREGROUND/ C BACKGROUND RELATIONSHIPS BETWEEN THE LINES, AND C THUS THE ORDER IN WHICH THEY SHOULD BE PLOTTED. C (XREF) AND (YREF) ARE THE COORDINATES, IN INCHES, C RELATIVE TO THE PLOTTER ORIGIN, TO BE USED AS THE ORIGIN C OF THE FIGURE. C (XLENTH) IS THE LENGTH, IN INCHES, TO WHICH THE PLOT IS C RESTRICTED. ANY POINT WHICH EXCEEDS THIS LIMIT, OR THE C LIMITS OF THE PAPER IN THE Y DIRECTION (NYPI), WILL BE C SET TO THAT LIMIT. C (MASK) IS AN INTEGER ARRAY OF 2*XLENTH*PIPI ENTRIES WHICH C IS USED TO STORE THE MASK. THE CONTENTS OF THIS ARRAY C SHOULD NOT BE ALTERED DURING THE PLOTING OF ANY GIVEN C FIGURE. C ALL PARAMETERS EXCEPT (MASK) AND (VERTEX) ARE RETURNED C UNCHANGED. C BETWEEN ANY TWO CALLS FOR THE SAME FIGURE, ANY PARAMETER C CAN BE MEANINGFULLY CHANGED EXCEPT (XLENTH), (MASK), AND C (VERTEX). INTEGER HIGH, OLDHI, OLDLOW DIMENSION XDATA(1), YDATA(1), ZDATA(1), MASK(1), * VERTEX(1) DATA INIT, JVXYZ, SPHI, STHETA/-1, -1, -1.0E99, * -1.0E99/ C INITIALIZATION PROCEDURES C INITIALIZATION PROCEDURE FOR A NEW FIGURE C TEST FOR SPECIAL MASK MODIFYING CALL IF (NLINE.EQ.0) GO TO 550 C DETERMINE IF INITIALIZATION IS REQUIRED IF (NLINE.NE.1) GO TO 20 C SET PLOTTER PARAMETERS PIPI = 100.0 NYPI = 1090 C RESET PLOTTER ORIGIN TO BOTTOM OF PLOT PAGE I = NYPI + 100 CALL IPLOT(0, -I, -3) C COMPUTE LENGTH OF PLOT PAGE IN INCREMENTS LIMITX = XLENTH*PIPI + 0.5 I = LIMITX + LIMITX C INITIALIZE MASKING ARRAY OVER THE LENGTH OF THE PLOT PAGE DO 10 K=1,I MASK(K) = INIT 10 CONTINUE INIT = -1 C SET THE NECESSARY INDICATORS FOR THE FIRST LINE OF A NEW C FIGURE INCI = -1 I = 0 C INPUT TYPE AND VERTEX INITIALIZATION C DETERMINE IF INITIALIZATION IS REQUIRED 20 IF (JVXYZ.EQ.IVXYZ) GO TO 70 C SET INDICATORS FOR TYPES OF INPUT DATA AND SAVING VERTICES JVXYZ = IVXYZ INDZ = 1 INDY = 1 INDX = 1 INDV = 1 IF (JVXYZ.LT.1000) GO TO 30 INDV = 2 JVXYZ = JVXYZ - 1000 30 IF (JVXYZ.LT.100) GO TO 40 INDX = 2 JVXYZ = JVXYZ - 100 40 IF (JVXYZ.LT.10) GO TO 50 INDY = 2 JVXYZ = JVXYZ - 10 50 IF (JVXYZ.LT.1) GO TO 60 INDZ = 2 60 JVXYZ = IVXYZ C ROTATION INITIALIZATION C DETERMINE IF INITIALIZATION IS REQUIRED 70 IF (PHI.EQ.SPHI .AND. THETA.EQ.STHETA) GO TO 80 C COMPUTE ROTATION FACTORS SPHI = SINF(0.0174532925*PHI) CPHI = COSF(0.0174532925*PHI) STHETA = SINF(0.0174532925*THETA) CTHETA = COSF(0.0174532925*THETA) A11 = CPHI A13 = -SPHI A21 = STHETA*SPHI A22 = CTHETA A23 = STHETA*CPHI SPHI = PHI STHETA = THETA C PROCESSING PROCEDURES C SET FLAG TO MOVE THROUGH THE DATA ARRAYS IN THE OPPOSITE C DIRECTION 80 INCI = -INCI C SET INDICATOR TO THE FIRST POINT TO BE PROCESSED IF (I.NE.0) I = NPNTS + 1 C LOOP TO PROCESS EACH POINT IN THE DATA ARRAYS DO 530 K=1,NPNTS C DATA CALCULATION I = I + INCI GO TO (90,100), INDX 90 X = XDATA + (I-1)*XSCALE GO TO 110 100 X = XDATA(I)*XSCALE 110 GO TO (120,130), INDY 120 Y = YDATA + (I-1)*YSCALE GO TO 140 130 Y = YDATA(I)*YSCALE 140 GO TO (150,160), INDZ 150 Z = ZDATA + (NLINE-1)*ZSCALE GO TO 170 160 Z = ZDATA(I)*ZSCALE C DATA ROTATION 170 XXX = A11*X + A13*Z + XREF XX = XXX IX = IROUND(XX*PIPI) YYY = A21*X + A23*Z + YREF YY = YYY + A22*Y IY = IROUND(YY*PIPI) C RESTRICT FIGURE TO PLOT PAGE IF (IX.LE.0) IX = 1 IF (IX.GT.LIMITX) IX = LIMITX IF (IY.LT.10) IY = 10 IF (IY.GT.NYPI) IY = NYPI IF (K.NE.1) GO TO 250 C (LOC) IS THE POSITION OF THE PREVIOUS POINT WITH RESPECT C TO THE MASK C +1 ABOVE THE MASK C 0 WITHIN THE LIMITS OF THE MASK C -1 BELOW THE MASK C PROCEDURE FOR INITIAL POINT OF EACH LINE C LOCATE INITIAL POINT WITH RESPECT TO THE MASK THEN C UPDATE THE MASK LOW = IX + IX HIGH = LOW - 1 MLOW = MASK(LOW) MHIGH = MASK(HIGH) IF (MHIGH-IY) 200, 210, 180 180 IF (MLOW-IY) 190, 230, 220 190 LOCOLD = 0 GO TO 240 200 MASK(HIGH) = IY IF (MLOW.EQ.-1) MASK(LOW) = IY 210 LOCOLD = +1 GO TO 240 220 MASK(LOW) = IY 230 LOCOLD = -1 C MOVE THE RAISED PEN TO THIS INITIAL POINT 240 CALL IPLOT(IX, IY, 3) JX = IX JY = IY IYREF = IY C STORE VERTICES IF REQUESTED IF (INDV.EQ.1) GO TO 530 INDEX = INCI + 6 VERTEX(INDEX) = XX VERTEX(INDEX+1) = YY VERTEX(INDEX+8) = XXX VERTEX(INDEX+9) = YYY IF (NLINE.NE.1) GO TO 530 VERTEX(1) = XX VERTEX(2) = YY VERTEX(9) = XXX VERTEX(10) = YYY GO TO 530 C SPECIAL CASE WHERE CHANGE IN X COORDINATE IS ZERO C A SPECIAL PROVISION IS MADE AT THIS POINT SO THAT A LINE C WILL NOT MASK ITSELF AS LONG AS THE X COORDINATE REMAINS C CONSTANT 250 IF (IX.NE.JX) GO TO 260 JY = IY GO TO 280 C COMPUTE CONSTANTS FOR LINEAR INTERPOLATION 260 YINC = FLOAT(IY-JY)/ABS(FLOAT(IX-JX)) INCX = (IX-JX)/IABS(IX-JX) YJ = JY C PREFORM LINEAR INTERPOLATION AT EACH INCREMENTAL STEP ON C THE X AXIS 270 JX = JX + INCX YJ = YJ + YINC JY = IROUND(YJ) C LOCATE THE CURRENT POINT WITH RESPECT TO THE MASK AT THAT C POINT THEN PLOT THE INCREMENT AS A FUNCTION OF THE C LOCATION OF THE PREVIOUS POINT WITH RESPECT TO ITS MASK LOW = JX + JX HIGH = LOW - 1 MLOW = MASK(LOW) MHIGH = MASK(HIGH) 280 IF (MHIGH-JY) 300, 300, 290 290 IF (MLOW-JY) 310, 320, 320 C THE CURRENT POINT IS ABOVE THE MASK 300 LOC = +1 IF (LOCOLD) 360, 370, 430 C THE CURRENT POINT IS WITHIN THE MASK 310 LOC = 0 IF (LOCOLD) 340, 350, 330 C THE CURRENT POINT IS BELOW THE MASK 320 LOC = -1 IF (LOCOLD) 510, 450, 440 C PLOT FROM ABOVE THE MASK TO WITHIN THE MASK 330 IF (MHIGH.LE.IYREF) CALL IPLOT(JX, MHIGH, 2) GO TO 350 C PLOT FROM BELOW THE MASK TO WITHIN THE MASK 340 IF (MLOW.GE.IYREF) CALL IPLOT(JX, MLOW, 2) C PLOT FROM WITHIN THE MASK TO WITHIN THE MASK 350 CALL IPLOT(JX, JY, 3) GO TO 520 C PLOT FROM BELOW THE MASK TO ABOVE THE MASK 360 IF (MLOW-IYREF) 370, 380, 380 C PLOT FROM WITHIN THE MASK TO ABOVE THE MASK 370 IF (MHIGH-IYREF) 400, 390, 390 380 CALL IPLOT(JX, MLOW, 2) 390 CALL IPLOT(JX, MHIGH, 3) GO TO 430 400 IF (MHIGH.EQ.-1) GO TO 430 OLDHI = HIGH - 2*INCX IF (MASK(OLDHI)-JY) 420, 420, 410 410 CALL IPLOT(JX, JY, 3) GO TO 430 420 CALL IPLOT(JX-INCX, MASK(OLDHI), 3) C PLOT FROM ABOVE THE MASK TO ABOVE THE MASK 430 MASK(HIGH) = JY IF (MLOW.EQ.-1) MASK(LOW) = JY CALL IPLOT(JX, JY, 2) GO TO 520 C PLOT FROM ABOVE THE MASK TO BELOW THE MASK 440 IF (MHIGH-IYREF) 460, 460, 450 C PLOT FROM WITHIN THE MASK TO BELOW THE MASK 450 IF (MLOW-IYREF) 470, 470, 480 460 CALL IPLOT(JX, MHIGH, 2) 470 CALL IPLOT(JX, MLOW, 3) GO TO 510 480 OLDLOW = LOW - 2*INCX IF (MASK(OLDLOW)-JY) 490, 500, 500 490 CALL IPLOT(JX, JY, 3) GO TO 510 500 CALL IPLOT(JX-INCX, MASK(OLDLOW), 3) C PLOT FROM BELOW THE MASK TO BELOW THE MASK 510 MASK(LOW) = JY CALL IPLOT(JX, JY, 2) 520 IYREF = JY LOCOLD = LOC IF (JX.NE.IX) GO TO 270 530 CONTINUE C RAISE PEN CALL IPLOT(JX, JY, 3) C STORE VERTICES IF REQUESTED IF (INDV.EQ.1) GO TO 540 INDEX = -INCI + 6 VERTEX(INDEX) = XX VERTEX(INDEX+1) = YY VERTEX(INDEX+8) = XXX VERTEX(INDEX+9) = YYY IF (NLINE.NE.1) GO TO 540 VERTEX(3) = XX VERTEX(4) = YY VERTEX(11) = XXX VERTEX(12) = YYY 540 I = I - 1 C RETURN TO CALLING PROGRAM RETURN C OPTION TO MODIFY THE MASKING TECHNIQUE TO BE USED ON THE C FOLLOWING FIGURE SO AS TO PLOT ONLY ABOVE ALL PREVIOUS C LINES. 550 INIT = 0 RETURN END SUBROUTINE FRAMER(IHCOR, VERTEX, MASK) C 10 C ROUTINE TO PLOT A FRAME ON THE PROJECTION OF A C 3-DIMENSIONAL FIGURE AS DRAWN BY PLOT3D. C INPUT PARAMETERS - C IHCOR NUMBER OF THE VERTEX OF THE FIGURE WHICH C APPEARS TO BE FURTHEST IN THE BACKGROUND C (MINUS Z DIRECTION). C VERTEX ARRAY CONTAINING THE COORDINATES OF THE C VERTICES OF THIS FIGURE AS RETURNED FROM C PLOT3D ON THE LAST CALL. C MASK ARRAY CONTAINING THE MASK FOR THIS FIGURE C AS RETURNED BY PLOT3D ON THE LAST CALL. C THE VERTICES OF THE FRAME ARE NUMBERED (1-4) IN THE SAME C ORDER AS THEIR COORDINATES APPEAR IN VERTEX. C THE MASK ARRAY IS ALTERED BY THIS ROUTINE, C BUT THE PLOTTER ORIGIN IS NOT MOVED. DIMENSION VERTEX(1), MASK(1), ARRAY(14) I = 2*IHCOR IF (I.LT.2) I = 2 IF (I.GT.8) I = 8 C THE VERTICES WHICH MAY BE HIDDEN C ARE DRAWN BY A CALL TO PLOT3D. ARRAY(1) = VERTEX(I-1) ARRAY(8) = VERTEX(I) ARRAY(2) = VERTEX(I+7) ARRAY(9) = VERTEX(I+8) ARRAY(4) = ARRAY(2) ARRAY(11) = ARRAY(9) ARRAY(6) = ARRAY(2) ARRAY(13) = ARRAY(9) ARRAY(7) = ARRAY(1) ARRAY(14) = ARRAY(8) I = I - 2 IF (I.EQ.0) I = 8 ARRAY(3) = VERTEX(I+7) ARRAY(10) = VERTEX(I+8) I = I + 4 IF (I.GT.8) I = I - 8 ARRAY(5) = VERTEX(I+7) ARRAY(12) = VERTEX(I+8) CALL PLOT3D(110, ARRAY, ARRAY(8), 0.0, 1.0, 1.0, 0.0, * 2, 7, 0.0, 0.0, 0.0, 0.0, 0.0, MASK, 0) C THE REMAINING VERTICES ARE DRAWN BY CALLS TO PLOT. CALL PLOT(VERTEX(I-1), VERTEX(I), 3) I = I - 2 DO 10 J=1,3 I = I + 2 IF (I.EQ.10) I = 2 CALL PLOT(VERTEX(I+7), VERTEX(I+8), 2) 10 CONTINUE CALL PLOT(VERTEX(I-1), VERTEX(I), 2) I = I - 2 IF (I.EQ.0) I = 8 CALL PLOT(VERTEX(I-1), VERTEX(I), 3) CALL PLOT(VERTEX(I+7), VERTEX(I+8), 2) RETURN END