home *** CD-ROM | disk | FTP | other *** search
/ Simtel MSDOS 1993 May / SIMTEL_0593.ISO / msdos / plot / 3d_plot.for < prev    next >
Text File  |  1988-01-18  |  9KB  |  269 lines
  1. Date: Tue, 29 Dec 87 13:50:53 EST
  2. From: Larry Granroth 319/335-1960 <"IOWASP::GRANROTH"@nssdca.GSFC.NASA.GOV>
  3. Subject: fortran code for 3d surface plotting
  4.  
  5.     Info-IBMPC Digest Volume 6, Issue 72, contained a request from
  6. <SIMPSON_P%MERCURY.ceo.dg.com@adam.dg.com> for a public domain 3-D
  7. plotting program.  I wrote a quick-and-dirty program which uses
  8. Tektronix PLOT10/TCS calls and ran on a VAX several years ago.  The
  9. algorithm effectively plots "slices" through a surface, so it isn't
  10. perfect, but it is faster than more sophisticated "web" algorithms.
  11. It should be fairly easy to convert to any FORTRAN system and
  12. substitute different plot drivers, so if anyone is interested, here is
  13. the source code:
  14.  
  15.       PROGRAM TEK3D
  16. C----------------------------------------------------------------------
  17. C               
  18. C       DEMONSTRATION PROGRAM TO DRIVE TEK 4014 TYPE TERMINAL
  19. C       MUST BE LINKED WITH PLOT10/TCS TYPE LIBRARY.
  20. C       DRAWS PERSPECTIVE SURFACE OF SINC(R) ON 40X40 GRID
  21. C       WITH R=SQRT(X*X+Y*Y)/2.
  22. C               
  23. C----------------------------------------------------------------------
  24.       DIMENSION Z(40,40)
  25.       COMMON /BLK3D/ IXCTR,IZCTR,IWIDE,IDEEP,IHIGH,VDIST
  26.       IXCTR=512 
  27.       IZCTR=390 
  28.       IWIDE=512 
  29.       IDEEP=512 
  30.       IHIGH=480 
  31.       VDIST=2048.
  32. C               
  33. H 67  528  528  0     0   528     0   1 1     0H        DO 20 J=1,40
  34.         Y=FLOAT(J)-20.5
  35.         YY=Y*Y
  36.         DO 10 I=1,40
  37.         X=FLOAT(I)-20.5
  38.         R=SQRT(X*X+YY)/2.
  39.         Z(I,J)=SIN(R)/R
  40.    10   CONTINUE
  41.    20   CONTINUE
  42. C               
  43.    99 TYPE '('' ENTER AZ,EL, AND IFLAG: '',$)'
  44.       ACCEPT *,AZ,EL,IFLAG
  45.       IF (IFLAG.LT.0) GOTO 9999
  46.       CALL PLT3D(Z,40,40,-.8,.8,AZ,EL,IFLAG)
  47.       PAUSE     
  48.       GOTO 99   
  49. C               
  50.  9999 STOP      
  51.       END       
  52. C----------------------------------------------------------------------
  53. C               
  54. C     RULED SURFACE VERSION OF PLT3D WRITTEN BY LARRY GRANROTH
  55. C               
  56. C     TEKTRONIX 4014 VERSION USES PLOT10/TCS ROUTINES
  57. C               
  58. C     VAX FORTRAN VERSION
  59. C               
  60. C     PARAMETERS:
  61. C     DATA     2-DIM REAL DATA ARRAY
  62. C     NX,NY    X AND Y INTEGER DIMENSIONS OF DATA ARRAY
  63. C     DMIN     DATA VALUE TO BE SCALED TO BOTTOM OF 3-D PLOT CUBE
  64. C     DMAX     VALUE TO BE SCALED TO TOP OF PLOT CUBE
  65. C     AZ       ANGLE (DEG) TO ROTATE PLOT CUBE CCW
  66. C     EL       ANGLE TO TILT PLOT CUBE TOWARD VIEWER
  67. C     IFLAG    1=X-SCAN ONLY, 2=Y-SCAN ONLY, OTHER=BOTH
  68. C               
  69. C     COMMON BLOCK:   /BLK3D/
  70. C     IXCTR    HORIZONTAL PIXEL COORDINATE TO CENTER PLOT CUBE
  71. C     IZCTR    VERTICLE COORDINATE  (16 BIT INTEGERS!)
  72. C     IWIDE    PLOT CUBE WIDTH IN PIXELS
  73. C     IDEEP    PLOT CUBE DEPTH
  74. C     IHIGH    PLOT CUBE HIGHT
  75. C     VDIST    VIEWING DISTANCE IN PIXELS  (FLOATING POINT)
  76. C               
  77. C     MODIFICATIONS:
  78. C     1-28-82  X OR Y SCAN OPTION AND TOTAL 360 AZ  -LJG
  79. C               
  80. C     8-26-82  PERSPECTIVE FORSHORTENING ADDED
  81. C              VDIST ADDED TO COMMON BLOCK          -LJG
  82. C
  83. C       This source code is placed in the public domain by
  84. C       Larry Granroth as of 29 December 1987.  I claim no
  85. C       responsibility for any problems associated with this
  86. C       code.
  87. C
  88. C----------------------------------------------------------------------
  89.       SUBROUTINE PLT3D(DATA,NX,NY,DMIN,DMAX,AZ,EL,IFLAG)
  90.       INTEGER MAXHID(1024), MINHID(1024)
  91.       LOGICAL PEN,NOHID,ABOV,BELO
  92.       DIMENSION DATA(NX,NY)
  93.       COMMON /BLK3D/ IXCTR,IZCTR,IWIDE,IDEEP,IHIGH,VDIST
  94.       COMMON /ROTRBK/ SINAZ,COSAZ,SINEL,COSEL
  95. C               
  96.       CALL INITT (1200)
  97. C               
  98.       XSCL=FLOAT(IWIDE)/FLOAT(NX)
  99. H 67  528  528  0     0   528     0   1 1     0H      YSCL=FLOAT(IDEEP)/FLOAT(NY)
  100.       IF (DMIN.EQ.DMAX) DMAX=DMIN+1.
  101.       ZSCL=FLOAT(IHIGH)/(DMAX-DMIN)
  102.       AZ=AZ*.01745329
  103.       EL=EL*.01745329
  104.       COSAZ=COS(AZ)
  105.       SINAZ=SIN(AZ)
  106.       COSEL=COS(EL)
  107.       SINEL=SIN(EL)
  108.       X0=(FLOAT(NX)-1.)/2.*XSCL
  109.       Y0=(FLOAT(NY)-1.)/2.*YSCL
  110.       D0=-(DMIN+(DMAX-DMIN)/2.)
  111.       IF (IFLAG.EQ.2) GOTO 5
  112.          IF (COSAZ.LT.0.) GOTO 3
  113.          LOY1=1 
  114.          LOY2=NY
  115.          IEX=1  
  116.          LIX1=2 
  117.          LIX2=NX
  118.          YST1=-Y0
  119.          XSC1=-X0
  120.          YSTI=YSCL
  121.          XSCI=XSCL
  122.          GO TO 4
  123.     3    LOY1=-NY
  124.          LOY2=-1
  125.          IEX=NX 
  126.          LIX1=1-NX
  127.          LIX2=-1
  128.          YST1=Y0
  129.          XSC1=X0
  130.          YSTI=-YSCL
  131.          XSCI=-XSCL
  132.     4    IF (IFLAG.EQ.1) GOTO 9
  133.     5    IF (SINAZ.LT.0.) GOTO 6
  134.          LOX1=1 
  135.          LOX2=NX
  136.          IEY=NY 
  137.          LIY1=1-NY
  138.          LIY2=-1
  139.          XST1=-X0
  140.          YSC1=Y0
  141.          XSTI=XSCL
  142.          YSCI=-YSCL
  143.          GO TO 8
  144.     6    LOX1=-NX
  145.          LOX2=-1
  146.          IEY=1  
  147.          LIY1=2 
  148.          LIY2=NY
  149.          XST1=X0
  150.          YSC1=-Y0
  151.          XSTI=-XSCL
  152.          YSCI=YSCL
  153.     8    IF (IFLAG.EQ.2) GOTO 100
  154.     9 DO 10 I=1,1024
  155.       MINHID(I)=959
  156.    10 MAXHID(I)=-960
  157.       Y1=YST1-YSTI
  158.       DO 50 LJ=LOY1,LOY2
  159.          J=IABS(LJ)
  160.          X1=XSC1
  161.          Y1=Y1+YSTI
  162.          Z1=(DATA(IEX,J)+D0)*ZSCL
  163.          CALL ROTR(X1,Y1,Z1,IX2,IZ2)
  164.          PEN=.FALSE.
  165. H 67  528  528  0     0   528     0   1 1     0H         IF (IZ2.GT.MAXHID(IX2).OR.IZ2.LT.MINHID(IX2)) PEN=.TRUE.
  166.          CALL MOVABS (IX2,IZ2)  ! MOVE PEN TO BEGINNING OF SCAN LINE
  167.          DO 40 LI=LIX1,LIX2
  168.          I=IABS(LI)
  169.          X1=X1+XSCI
  170.          IX1=IX2
  171.          IZ1=IZ2
  172.          Z1=(DATA(I,J)+D0)*ZSCL
  173.          CALL ROTR(X1,Y1,Z1,IX2,IZ2)
  174.          NOHID=.TRUE.
  175.          IF (IX2.EQ.IX1) GOTO 40
  176.          SLOPE=FLOAT(IZ2-IZ1)/FLOAT(IX2-IX1)
  177.          XX=0.  
  178.          IZ=IZ1 
  179.          IX=IX1+1
  180.          DO 30 K=IX,IX2
  181.             XX=XX+1.
  182.             LAST=IZ1
  183.             IZ1=INT(SLOPE*XX+.5)+IZ
  184.             ABOV=(IZ1.GT.MAXHID(K))
  185.             BELO=(IZ1.LT.MINHID(K))
  186.             IF (ABOV.OR.BELO) GOTO 25
  187.             NOHID=.FALSE.
  188.             IF (PEN) CALL DRWABS(K-1,LAST)  ! DRAW TO LAST UNHIDDEN POINT
  189.             PEN=.FALSE.
  190.             GOTO 30
  191.    25       IF (ABOV) MAXHID(K)=IZ1
  192.             IF (BELO) MINHID(K)=IZ1
  193.             IF (PEN) GOTO 30
  194.             CALL MOVABS(K,IZ1)   ! MOVE PEN
  195.             PEN=.TRUE.
  196.    30       CONTINUE
  197.    40    IF (NOHID) CALL DRWABS(IX2,IZ2)
  198.    50    CONTINUE
  199.       IF (IFLAG.EQ.1) GOTO 150
  200.   100 DO 110 I=1,1024
  201.       MINHID(I)=959
  202.   110 MAXHID(I)=-960
  203.       X1=XST1-XSTI
  204.       DO 150 LJ=LOX1,LOX2
  205.          J=IABS(LJ)
  206.          Y1=YSC1
  207.          X1=X1+XSTI
  208.          Z1=(DATA(J,IEY)+D0)*ZSCL
  209.          CALL ROTR(X1,Y1,Z1,IX2,IZ2)
  210.          PEN=.FALSE.
  211.          IF (IZ2.GT.MAXHID(IX2).OR.IZ2.LT.MINHID(IX2)) PEN=.TRUE.
  212.          CALL MOVABS(IX2,IZ2)   ! MOVE PEN TO BEGINNING OF SCAN LINE
  213.          DO 140 LI=LIY1,LIY2
  214.          I=IABS(LI)
  215.          Y1=Y1+YSCI
  216.          IX1=IX2
  217.          IZ1=IZ2
  218.          Z1=(DATA(J,I)+D0)*ZSCL
  219.          CALL ROTR(X1,Y1,Z1,IX2,IZ2)
  220.          NOHID=.TRUE.
  221.          IF (IX2.EQ.IX1) GOTO 140
  222.          SLOPE=FLOAT(IZ2-IZ1)/FLOAT(IX2-IX1)
  223.          XX=0.  
  224.          IZ=IZ1 
  225.          IX=IX1+1
  226.          DO 130 K=IX,IX2
  227.             XX=XX+1.
  228.             LAST=IZ1
  229.             IZ1=INT(SLOPE*XX+.5)+IZ
  230.             ABOV=(IZ1.GT.MAXHID(K))
  231. H 67  528  528  0     0   528     0   1 1     0H            BELO=(IZ1.LT.MINHID(K))
  232.             IF (ABOV.OR.BELO) GOTO 125
  233.             NOHID=.FALSE.
  234.             IF (PEN) CALL DRWABS(K-1,LAST)   ! DRAW TO LAST UNHIDDEN POINT
  235.             PEN=.FALSE.
  236.             GOTO 130
  237.   125       IF (ABOV) MAXHID(K)=IZ1
  238.             IF (BELO) MINHID(K)=IZ1
  239.             IF (PEN) GOTO 130
  240.             CALL MOVABS(K,IZ1)   ! MOVE PEN
  241.             PEN=.TRUE.
  242.   130       CONTINUE
  243.   140    IF (NOHID) CALL DRWABS(IX2,IZ2)
  244.   150    CONTINUE
  245.       CALL FINITT(0,780)
  246.       RETURN    
  247.       END       
  248. C----------------------------------------------------------------------
  249.       SUBROUTINE ROTR(X1,Y1,Z1,IP,KP)
  250.       COMMON /BLK3D/ IXCTR,IZCTR,IWIDE,IDEEP,IHIGH,VDIST
  251.       COMMON /ROTRBK/ SINAZ,COSAZ,SINEL,COSEL
  252.       XROT =X1*COSAZ-Y1*SINAZ
  253.       YROT =X1*SINAZ+Y1*COSAZ
  254.       ZTILT=YROT*SINEL+Z1*COSEL
  255.       PERSP=VDIST/(VDIST+(YROT*COSEL-Z1*SINEL))
  256.       IP   =INT(XROT *PERSP+0.5)+IXCTR
  257.       KP   =INT(ZTILT*PERSP+0.5)+IZCTR
  258.       RETURN    
  259.       END       
  260.  
  261. Various ways to reach me include:
  262. SPAN      IOWASP::GRANROTH
  263. ARPAnet   GRANROTH%IOWASP.SPAN@NSSDC.GSFC.NASA.GOV
  264.           or                  @STAR.STANFORD.EDU
  265.           or                  @VLSI.JPL.NASA.GOV
  266.           or                  @SDS
  267. BITNET    CCOLJGPG@UIAMVS
  268. TELEMAIL  LGRANROTH/J.P.L.
  269.