home *** CD-ROM | disk | FTP | other *** search
/ RISC DISC 2 / RISC_DISC_2.iso / pd_share / program / fortran77_210 / acmtoms / f77 / surface475 < prev   
Encoding:
Text File  |  1992-01-04  |  19.2 KB  |  510 lines
  1. ***********************  COPYRIGHT NOTICE  ********************************
  2. *  The material in this library is copyrighted by the ACM, which grants   *
  3. *  general permission to distribute provided the copies are not made for  *
  4. *  direct commercial advantage.  For details of the copyright and         *
  5. *  dissemination agreement, consult the current issue of TOMS.            *
  6. ***************************************************************************
  7. C
  8. C          *** from netlib, Thu Jan  2 14:42:23 GMT 1992 ***
  9. C     ALGORITHM 475 COLLECTED ALGORITHMS FROM ACM.
  10. C     ALGORITHM APPEARED IN COMM. ACM, VOL. 17, NO. 03,
  11. C     P. 152.
  12. C             Archimedes graphics added KM Crennell 2 Jan 92
  13. C
  14. C  ********           Needs PD-F77 graphics library       ********
  15. C
  16.       PROGRAM ACMTEST                               
  17. C DEMONSTRATION PROGRAM                             
  18.       DIMENSION EYE(3),S(4),ST1(80,80,2),IS2(3,160)
  19.       DIMENSION IOBJ(80,80)                        
  20.       CHARACTER *22 SCRFIL
  21. C                    unit number for scratch file
  22.       DATA IUS/9/
  23. C                     name of scratch file
  24.       DATA SCRFIL/'SCRATCH'/
  25. C                       initialise Graphics   KM Crennell
  26.       CALL GRINIT(XMM,YMM)
  27. C           Xmm,Ymm is size of graphics screen, force a square area
  28.       XSQ=XMM
  29.       IF(XMM.GT.YMM)XSQ=YMM
  30. C USE WHOLE FRAME    
  31.       S(1)=0.                                      
  32.       S(2)=XSQ
  33.       S(3)=0.                                      
  34.       S(4)=XSQ                                   
  35. C SET EYE POSITION                                 
  36.       EYE(1)=250.                                  
  37.       EYE(2)=150.                                  
  38.       EYE(3)=100.                                  
  39. C                        set up scratch file
  40.       OPEN(UNIT=IUS,FILE=SCRFIL,ACCESS='SEQUENTIAL',err=80,
  41.      1    FORM='UNFORMATTED',STATUS='OLD')
  42.       ICON=4
  43.       IF(ICON.EQ.4)GOTO 202
  44.       GOTO 200
  45.   80  WRITE(*,*)' error opening scratch file ',SCRFIL
  46.       GOTO60
  47. C                        INITIALIZE PACKAGE
  48.   200 CALL INIT3D(EYE,80,80,80,ST1,3,160,IS2,IUS,S)
  49. C                        CREATE AND PLOT 1st TEST OBJECT 
  50.   202 CALL TEST1(80,80,ST1,3,160,160,IS2,IUS,S,IOBJ,80)
  51. C                        ADVANCE TO THE NEXT FRAME.                           
  52. C  204 CALL GEMPTY
  53.       WRITE(*,*)' 1st pic done, type any number to continue'
  54. C                          pause and look at picture
  55.       READ(*,*)KMC   
  56.       CALL CLS
  57. C                        clear the screen
  58. C      CALL GCLEAR
  59. C A SECOND PICTURE WILL NOW BE CALLED USING THE SAME SIZE   
  60. C ARRAYS AND EYE POSITION.  THIS MEANS THE CALL TO INIT3D,  
  61. C THE BIGGEST TIME CONSUMER, CAN BE SKIPPED IF THE FOLLOWING
  62. C FOUR LINES ARE INCLUDED.                                  
  63.    25 REWIND IUS                                              
  64.       DO 5 I=1,3                                            
  65.         DO 5 J=1,160                                        
  66.     5     IS2(I,J)=0  
  67. C                     now make and plot 2nd picture
  68. C                    set colour depending on the level
  69.       IF(ICON.EQ.4)THEN
  70.         KOL=1
  71.       ELSE
  72.         KOL=ICON
  73.       ENDIF
  74. C      CALL GCOL(KOL)
  75.       CALL TEST2(80,80,ST1,3,160,160,IS2,IUS,S,IOBJ,80,ICON)
  76. C FLUSH PLOT BUFFER                                         
  77. C      CALL FRAME 
  78.       WRITE(*,*)' end of 2nd Pic ICON=',ICON,' enter new value '
  79.       READ(*,*)KMC
  80.       IF(KMC.LE.0)GOTO50
  81.       ICON=KMC
  82.       WRITE(*,*)' enter 0 to leave previous plot, 1 to clear it'
  83.       READ(*,*)KMC
  84.       IF(KMC.EQ.1)CALL CLS
  85.       GOTO 25
  86. C                        close the scratch file Should delete really
  87.    50 CLOSE(IUS)
  88. C                         stop the Graphics         
  89.    60 CALL GREND
  90.       STOP                                                  
  91.       END       
  92.       SUBROUTINE TEST1(NV,NW,ST1,LX,NX,NY,IS2,IUS,S,IOBJ,MV)
  93. C                            CREATE AND PLOT 1st TEST OBJECT  
  94.       DIMENSION ST1(NV,NW,2),IS2(LX,NY),IOBJ(NV,NW)
  95.       DO 4 I=1,80                                  
  96.         A=(I-50)**2                                
  97.         DO 3 J=1,80                                
  98.           C=(J-25)**2                              
  99.           D=IABS(J-63)+IABS(I-25)                  
  100.           DO 3 K=1,80                              
  101. C FLOOR                                            
  102.             IF(K.EQ.1) GO TO 1                     
  103. C BALL                                             
  104.             IF(SQRT(A+C+(FLOAT(K)-25.)**2).LE.25.) GO TO 1
  105. C POINT
  106.             IF(D.GT.FLOAT(80-K)*.1875) GO TO 2
  107. C                    object fills this voxel
  108.     1       IOBJ(J,K)=1                       
  109.             GO TO 3                         
  110. C                   no object in this voxel  
  111.     2       IOBJ(J,K)=0                       
  112.     3       CONTINUE 
  113.         CALL DANDR(NV,NW,ST1,LX,NX,NY,IS2,IUS,S,IOBJ,MV,IER)
  114. C                             allow for eof on scratch file
  115.         IF(IER.NE.0)RETURN
  116.     4   CONTINUE
  117.       RETURN
  118.       END
  119.       SUBROUTINE TEST2(NV,NW,ST1,LX,NX,NY,IS2,IUS,S,IOBJ,MV,ICON)
  120. C                            CREATE AND PLOT 2nd TEST OBJECT  
  121.       DIMENSION ST1(NV,NW,2),IS2(LX,NY),IOBJ(NV,NW)
  122. C                 THIS PICTURE WILL BE THE T=ICON CONTOUR SURFACE OF           
  123. C           T=1/SQRT(U*U+V*V+W*W)+(.5-V)**2/SQRT(U*U+V*V).            
  124.       T=ICON
  125.       DO 9 I=1,80                                           
  126.         U=(40.5-FLOAT(I))/79.                               
  127.         UU=U*U                                              
  128.         DO 8 J=1,80                                         
  129.           V=(FLOAT(J)-40.5)/79.                             
  130.           VV=V*V                                            
  131.           A=1./SQRT(UU+VV)                                  
  132.           DO 8 K=1,80                                       
  133. C           THE FOLLOWING CARD ADDS AXES.                             
  134.             IF(I*J.EQ.1.OR.I*K.EQ.1.OR.J*K.EQ.1) GO TO 6    
  135.             W=(FLOAT(K)-40.5)/79.
  136. C                                                *  contour level 
  137.             IF(1./SQRT(UU+VV+W*W)+(.5-V)**2*A.LE.T) GO TO 7
  138.     6       IOBJ(J,K)=1                                     
  139.             GO TO 8                                         
  140.     7       IOBJ(J,K)=0                                     
  141.     8     CONTINUE                     
  142.           CALL DANDR(NV,NW,ST1,LX,NX,NY,IS2,IUS,S,IOBJ,MV,IER)
  143.         IF(IER.NE.0)RETURN
  144.     9 CONTINUE
  145.       RETURN
  146.       END
  147.       SUBROUTINE GRINIT(XMM,YMM)
  148. C                         graphics initialisation 
  149.       CALL MODE(27)
  150.       XMM=1280
  151.       YMM=960
  152.       RETURN
  153.       END
  154.       SUBROUTINE GREND
  155. C                         graphics close
  156.       RETURN
  157.       END
  158.       SUBROUTINE  CLINE(X1,Y1,X2,Y2)
  159. C              IS ASSUMED TO DRAW A LINE  FROM (X1,Y1) TO (X2,Y2)
  160. C            what co-ordinates? assume plotter ones. KM Crennell
  161.       IX1=X1
  162.       IY1=Y1
  163.       IX2=X2
  164.       IY2=Y2
  165.       CALL LINE(IX1,IY1,IX2,IY2)
  166.       RETURN
  167.       END
  168.       SUBROUTINE INIT3D(EYE,NU,NV,NW,ST1,LX,NY,IS2,IU,S)    
  169.       DIMENSION EYE(3),ST1(NV,NW,2),IS2(LX,NY),S(4)
  170. C
  171. C BY THOMAS WRIGHT
  172. C COMPUTING FACILITY
  173. C THE NATIONAL CENTER FOR ATMOSPHERIC RESEARCH
  174. C BOULDER, COLORADO 80302
  175. C NCAR IS SPONSORED BY THE NATIONAL SCIENCE FOUNDATION.
  176. C
  177. C THE METHOD IS DESCRIBED IN DETAIL IN - A ONE-PASS HIDDEN-
  178. C LINE REMOVER FOR COMPUTER DRAWN THREE-SPACE OBJECTS. PROC
  179. C 1972 SUMMER COMPUTER SIMULATION CONFERENCE, 261-267, 1972.
  180. C
  181. C THIS VERSION IS FOR USE ON CDC 6000 OR 7000 COMPUTERS.
  182. C
  183. C THIS PACKAGE OF ROUTINES PLOTS 3-DIMENSIONAL OBJECTS WITH
  184. C HIDDEN PARTS NOT SHOWN.  OBJECTS ARE STORED IN AN ARRAY,
  185. C WITH THE POSITION IN THE ARRAY CORRESPONDING TO A LOCATION
  186. C IN 3-SPACE AND THE VALUE OF THE ARRAY ELEMENT TELLING IF
  187. C ANY OBJECT IS PRESENT AT THE LOCATION.
  188. C
  189. C INIT3D IS AN INITIALIZATION ROUTINE FOR THIS PACKAGE.  IT
  190. C IS CALLED, THEN A SEQUENCE OF CALLS ARE MADE TO DANDR TO
  191. C PRODUCE A PICTURE.
  192. C EYE   AN ARRAY 3 LONG CONTAINING THE U, V, AND W COORDI-
  193. C       NATES OF THE EYE POSITION.  OBJECTS ARE CONSIDERED
  194. C       TO BE IN A BOX WITH 2 EXTREME CORNERS AT (1,1,1) AND
  195. C       (NU,NV,NW).  THE EYE POSITION MUST HAVE POSITIVE
  196. C       COORDINATES AWAY FROM THE COORDINATE PLANES U=0,
  197. C       V=0, AND W=0.  WHILE GAINING EXPERIENCE WITH THE
  198. C       PACKAGE, USE EYE(1)=5*NU, EYE(2)=4*NV, EYE(3)=3*NW.
  199. C NU    U DIRECTION LENGTH OF THE BOX CONTAINING THE OBJECTS
  200. C NV    V DIRECTION LENGTH OF THE BOX CONTAINING THE OBJECTS
  201. C NW    W DIRECTION LENGTH OF THE BOX CONTAINING THE OBJECTS
  202. C ST1   A SCRATCH ARRAY AT LEAST NV*NW*2 WORDS LONG.
  203. C LX    FIRST DIMENSION OF A SCRATCH ARRAY, IS2, USED BY THE
  204. C       PACKAGE FOR REMEMBERING WHERE IT SHOULD NOT DRAW.
  205. C       LX=1+NX/NBPW.  SEE DANDR COMMENTS FOR NX AND NBPW.
  206. C NY    SECOND DIMENSION OF IS2.  SEE DANDR COMMENTS.
  207. C IS2   A SCRATCH ARRAY AT LEAST LX*NY WORDS LONG.
  208. C IU    UNIT NUMBER OF SCRATCH FILE FOR THE PACKAGE.  ST1
  209. C       WILL BE WRITTEN NU TIMES ON THIS FILE.
  210. C S     AN ARRAY 4 LONG WHICH CONTAINS THE COORDINATES OF
  211. C       THE AREA WHERE THE PICTURE IS TO BE DRAWN.  THAT IS,
  212. C       ALL PLOTTING COORDINATES GENERATED WILL BE BOUNDED
  213. C       AS FOLLOWS-- X COORDINATES WILL BE BETWEEN S(1) AND
  214. C       S(2), Y COORDINATES WILL BE BETWEEN S(3) AND S(4).
  215. C       TO PREVENT DISTORTION, HAVE S(2)-S(1)=S(4)-S(3).
  216. C
  217. C IF SEVERAL PICTURES ARE TO BE DRAWN WITH THE SAME SIZE
  218. C ARRAYS AND EYE POSITION AND THE USER REWINDS IU AND FILLS
  219. C IS2 WITH ZEROES, INIT3D NEED NOT BE CALLED FOR OTHER THAN
  220. C THE FIRST PICTURE.
  221. C
  222. C SET UP TRANSFORMATION ROUTINE FOR THIS LINE OF SIGHT.
  223.       U=NU
  224.       V=NV
  225.       W=NW
  226.       CALL SETORG(U*.5,V*.5,W*.5,EYE(1),EYE(2),EYE(3))
  227. C FIND EXTREMES IN TRANSFORMED SPACE.
  228.       CALL PERSPC(1.,1.,W,D,YT,D)
  229.       CALL PERSPC(U,V,1.,D,YB,D)
  230.       CALL PERSPC(U,1.,1.,XL,D,D)
  231.       CALL PERSPC(1.,V,1.,XR,D,D)
  232. C ADJUST EXTREMES TO PREVENT DISTORTION WHEN GOING FROM
  233. C TRANSFORMED SPACE TO PLOTTER SPACE.
  234.       DIF=(XR-XL-YT+YB)*.5
  235.       IF(DIF) 1,3,2
  236.     1 XL=XL+DIF
  237.       XR=XR-DIF
  238.       GO TO 3
  239.     2 YB=YB-DIF
  240.       YT=YT+DIF
  241.     3 REWIND IU                    
  242. C FIND THE PLOTTER COORDINATES OF THE 3-SPACE LATTICE POINTS
  243.       C1=.9*(S(2)-S(1))/(XR-XL)
  244.       C2=.05*(S(2)-S(1))+S(1)
  245.       C3=.9*(S(4)-S(3))/(YT-YB)
  246.       C4=.05*(S(4)-S(3))+S(3)     
  247.       DO 5 I=1,NU
  248.         U=NU+1-I
  249.         DO 4 J=1,NV
  250.           V=J
  251.           DO 4 K=1,NW
  252.             CALL PERSPC(U,V,FLOAT(K),X,Y,D)
  253.             ST1(J,K,1)=C1*(X-XL)+C2
  254.     4       ST1(J,K,2)=C3*(Y-YB)+C4
  255. C WRITE THEM ON UNIT IU.
  256.     5   WRITE(IU) ST1
  257.       REWIND IU
  258. C ZERO OUT ARRAY WHERE VISIBILITY IS REMEMBERED.
  259.       DO 6 J=1,NY
  260.         DO 6 I=1,LX
  261.     6     IS2(I,J)=0
  262.       RETURN
  263.       END
  264.       SUBROUTINE SETORG(X,Y,Z,XT,YT,ZT)           
  265. C
  266. C THIS ROUTINE IMPLEMENTS THE 3-SPACE TO 2-SPACE TRANSFOR-
  267. C MATION BY KUBER, SZABO AND GIULIERI, THE PERSPECTIVE
  268. C REPRESENTATION OF FUNCTIONS OF TWO VARIABLES. J. ACM 15,
  269. C 2, 193-204,1968.
  270. C SETORG ARGUMENTS
  271. C X,Y,Z    ARE THE 3-SPACE COORDINATES OF THE INTERSECTION
  272. C          OF THE LINE OF SIGHT AND THE IMAGE PLANE.  THIS
  273. C          POINT CAN BE THOUGHT OF AS THE POINT LOOKED AT.
  274. C XT,YT,ZT ARE THE 3-SPACE COORDINATES OF THE EYE POSITION.
  275. C
  276. C PERSPC ARGUMENTS
  277. C X,Y,Z    ARE THE 3-SPACE COORDINATES OF A POINT TO BE
  278. C          TRANSFORMED.
  279. C XT,YT    THE RESULTS OF THE 3-SPACE TO 2-SPACE TRANSFOR-
  280. C          MATION.
  281. C ZT       NOT USED.
  282. C
  283. C STORE THE PARAMETERS OF THE SETORG CALL FOR USE WHEN
  284. C PERSPC IS CALLED.
  285.       AX=X
  286.       AY=Y
  287.       AZ=Z
  288.       EX=XT
  289.       EY=YT
  290.       EZ=ZT
  291. C AS MUCH COMPUTATION AS POSSIBLE IS DONE DURING EXECUTION
  292. C OF SETORG SINCE PERSPC IS CALLED THOUSANDS OF TIMES FOR
  293. C EACH CALL TO SETORG.
  294.       DX=AX-EX
  295.       DY=AY-EY
  296.       DZ=AZ-EZ
  297.       D=SQRT(DX*DX+DY*DY+DZ*DZ)
  298.       COSAL=DX/D
  299.       COSBE=DY/D
  300.       COSGA=DZ/D
  301.       AL=ACOS(COSAL)
  302.       BE=ACOS(COSBE)
  303.       GA=ACOS(COSGA)
  304.       SINGA=SIN(GA)
  305. C THE 3-SPACE POINT LOOKED AT IS TRANSFORMED INTO (0,0) OF
  306. C THE 2-SPACE.  THE 3-SPACE Z AXIS IS TRANSFORMED INTO THE
  307. C 2-SPACE Y AXIS.  IF THE LINE OF SIGHT IS CLOSE TO PARALLEL
  308. C TO THE 3-SPACE Z AXIS, THE 3-SPACE Y AXIS IS CHOSEN (IN-
  309. C STEAD OF THE 3-SPACE Z AXIS) TO BE TRANSFORMED INTO THE
  310. C 2-SPACE Y AXIS.
  311.       IF(SINGA.LT.0.0001) GO TO 1
  312.       R=1./SINGA
  313.       ASSIGN 2 TO JUMP
  314.       RETURN
  315.     1 SINBE=SIN(BE)
  316.       R=1./SINBE
  317.       ASSIGN 3 TO JUMP
  318.       RETURN
  319. C********************  ENTRY PERSPC  ***********************
  320.       ENTRY PERSPC(X,Y,Z,XT,YT,ZT)
  321.       Q=D/((X-EX)*COSAL+(Y-EY)*COSBE+(Z-EZ)*COSGA)
  322.       GO TO JUMP,(2,3)
  323.     2 XT=((EX+Q*(X-EX)-AX)*COSBE-(EY+Q*(Y-EY)-AY)*COSAL)*R
  324.       YT=(EZ+Q*(Z-EZ)-AZ)*R
  325.       RETURN
  326.     3 XT=((EZ+Q*(Z-EZ)-AZ)*COSAL-(EX+Q*(X-EX)-AX)*COSGA)*R
  327.       YT=(EY+Q*(Y-EY)-AY)*R
  328.       RETURN
  329.       END
  330.       SUBROUTINE DANDR(NV,NW,ST1,LX,NX,NY,IS2,IU,S,IOBJS,MV,IER)
  331.       DIMENSION ST1(NV,NW,2),IS2(LX,NY),S(4),IOBJS(MV,NW)
  332. C
  333. C THIS ROUTINE IS CALLED NU TIMES, EACH CALL PROCESSING THE
  334. C PART OF THE PICTURE AT U=NU+1-I WHERE I IS THE NUMBER OF
  335. C THE CALL TO DANDR.  THAT IS, THE PART OF THE PICTURE AT
  336. C U=NU IS PROCESSED DURING THE FIRST CALL, THE PART OF THE
  337. C PICTURE AT U=NU-1 IS PROCESSED DURING THE SECOND CALL, AND
  338. C SO ON UNTIL THE PART OF THE PICTURE AT U=1 IS PROCESSED
  339. C DURING THE LAST CALL.
  340. C NV    SEE INIT3D COMMENTS.
  341. C NW    SEE INIT3D COMMENTS.
  342. C ST1   SEE INIT3D COMMENTS.
  343. C LX    THE NUMBER OF WORDS NEEDED TO HOLD NX BITS.  ALSO,
  344. C       THE FIRST DIMENSION OF IS2.
  345. C NX    NUMBER OF CELLS IN THE X DIRECTION OF A MODEL OF THE
  346. C       IMAGE PLANE.  A SILHOUETTE OF THE PARTS OF THE PIC-
  347. C       TURE PROCESSED SO FAR IS STORED IN THIS MODEL. LINES
  348. C       TO BE DRAWN ARE TESTED FOR VISIBILITY BY EXAMINING
  349. C       THE SILHOUETTE.  LINES IN THE SILHOUETTE ARE HIDDEN.
  350. C       LINES OUT OF THE SILHOUETTE ARE VISIBLE.  THE SOLU-
  351. C       TION IS APPROXIMATE BECAUSE THE SILHOUETTE IS NOT
  352. C       FORMED EXACTLY.  SEE IS2 COMMENT BELOW.
  353. C NY    NUMBER OF CELLS IN THE Y DIRECTION OF THE MODEL OF
  354. C       THE IMAGE PLANE.  ALSO THE SECOND DIMENSION OF IS2.
  355. C IS2   AN ARRAY TO HOLD THE IMAGE PLANE MODEL.  IT IS
  356. C       DIMENSIONED LX BY NY.  THE MODEL IS NX BY NY AND
  357. C       PACKED DENSELY.  IF HIDDEN LINES ARE DRAWN, DECREASE
  358. C       NX AND NY (AND LX IF POSSIBLE).  IF VISIBLE LINES
  359. C       ARE LEFT OUT OF THE PICTURE, INCREASE NX AND NY (AND
  360. C       LX IF NEED BE).  AS A GUIDE, SOME EXAMPLES SHOWING
  361. C       SUCCESSFUL CHOICES ARE LISTED
  362. C          GIVEN  NU  NV  NW   RESULTING NX  NY FROM TESTING
  363. C                100 100  60            200 200
  364. C                 60  60  60            110 110
  365. C                 40  40  40             75  75
  366. C IU    SEE INIT3D COMMENTS.
  367. C IOBJS A NV BY NW ARRAY (WITH ACTUAL FIRST DIMENSION MV IN
  368. C       THE CALLING PROGRAM) DESCRIBING THE OBJECT.  IF THIS
  369. C       IS CALL NUMBER I TO DANDR, THE PART OF THE PICTURE
  370. C       AT U=NU+1-I IS TO BE PROCESSED.  IOBJS DEFINES THE
  371. C       OBJECTS TO BE DRAWN IN THE FOLLOWING MANNER --
  372. C       IOBJS(J,K)=1 IF ANY OBJECT CONTAINS THE POINT
  373. C       (NU+1-I,J,K) AND IOBJS(J,K)=0 OTHERWISE.
  374. C MV    ACTUAL FIRST DIMENSION OF IOBJS IN THE CALLING
  375. C       PROGRAM.
  376. C
  377. C************** MACHINE DEPENDANT CONSTANTS ****************
  378. C NBPW NUMBER OF BITS PER WORD
  379. C MASK AN ARRAY NBPW LONG.  MASK(I)=2**(I-1), I=1,2,...,NBPW
  380. C CDC 6000 OR 7000 VERSION changed to 32 bits for VAX  KM Crennell
  381.       DIMENSION MASK(32)
  382.       DATA NBPW/32/
  383.       DATA MASK/?I01,?I02,?I04,?I08,?I10,?I20,?I40,?I80,?I0100,
  384.      *?I0200,?I0400,?I0800,?I1000,?I2000,?I4000,?I8000,?I010000,
  385.      *  ?I020000,?I040000,?I080000,?I100000,?I200000,?I400000,
  386.      *  ?I800000,?I01000000,?I02000000,?I04000000,?I08000000,
  387.      *?I10000000,?I20000000,?I40000000,?I80000000/
  388. C
  389.       IER=0
  390.       ASSIGN 12 TO IRET
  391. C RX AND RY ARE USED TO MAP PLOTTER COORDINATES INTO THE
  392. C IMAGE PLANE MODEL.
  393.       RX=(FLOAT(NX)-1.)/(S(2)-S(1))
  394.       RY=(FLOAT(NY)-1.)/(S(4)-S(3))
  395. C READ THE RELATIVE PLOTTER COORDINATES OF THE LATTICE
  396. C POINTS FROM UNIT IU.
  397.       READ(IU,END=100) ST1
  398. C DX, DY AND DZ ARE USED TO FIND REQUIRED COORDINATES OF
  399. C NON-LATTICE POINTS.
  400.       NVD2=NV/2
  401.       NWD2=NW/2
  402.       DX=(ST1(NV,NWD2,1)-ST1(1,NWD2,1))*.5/(FLOAT(NV)-1.)
  403.       DY=(ST1(1,NWD2,2)-ST1(NV,NWD2,2))*.5/(FLOAT(NV)-1.)
  404.       DZ=(ST1(NVD2,NW,2)-ST1(NVD2,1,2))*.5/(FLOAT(NW)-1.)
  405. C SLOPE IS USED TO DEFORM THE IMAGE PLANE MODEL SO THAT
  406. C LINES OF CONSTANT Y OF THE IMAGE MODEL HAVE THE SAME
  407. C SLOPE AS LINES OF CONSTANT U AND W IN THE PICTURE.  THIS
  408. C IMPROVES THE PICTURE.
  409.       SLOPE=DY/DX
  410. C THE FOLLOWING LOOPS THROUGH STATEMENT 12 GENERATE THE .5
  411. C CONTOUR LINES IN 2-SPACE FOR THE ARRAY IOBJS (WHICH CON-
  412. C TAINS ONLY ZEROES AND ONES), TESTS THE LINES FOR VISIBIL-
  413. C ITY, AND CALLS A ROUTINE TO PLOT THE VISIBLE LINES.
  414.       DO 12 I=2,NV
  415.         JUMP=IOBJS(I-1,1)*8+IOBJS(I,1)*4+1
  416.         DO 12 J=2,NW
  417.           X=ST1(I,J,1)
  418.           Y=ST1(I,J,2)
  419. C DECIDE WHICH OF THE 16 POSSIBILITIES THIS IS.
  420.           JUMP=(JUMP-1)/4+IOBJS(I-1,J)*8+IOBJS(I,J)*4+1
  421.           GO TO (12,2,4,5,7,8,3,10,10,1,8,7,5,4,2,12),JUMP
  422. C GOING TO 1 MEANS JUMP=10 WHICH MEANS ONLY THE LOWER-RIGHT
  423. C AND UPPER-LEFT ELEMENTS OF THIS CELL ARE SET TO 1.
  424. C TWO LINES SHOULD BE DRAWN, A DIAGONAL CONNECTING THE
  425. C MIDDLE OF THE BOTTOM TO THE MIDDLE OF THE RIGHT SIDE OF
  426. C THE CELL (LOWER-RIGHT LINE), AND A DIAGONAL CONNECTING THE
  427. C MIDDLE OF THE LEFT SIDE TO THE MIDDLE OF THE TOP (UPPER-
  428. C LEFT LINE) OF THE CELL.
  429.     1     ASSIGN 9 TO IRET
  430. C LOWER-RIGHT LINE
  431.     2     X1=X
  432.           Y1=Y-DZ
  433.           X2=X+DX
  434.           Y2=Y-DY
  435.           GO TO 11
  436. C LOWER-LEFT AND UPPER-RIGHT
  437.     3     ASSIGN 6 TO IRET
  438. C LOWER-LEFT
  439.     4     X1=X
  440.           Y1=Y-DZ
  441.           X2=X-DX
  442.           Y2=Y+DY
  443.           GO TO 11
  444. C HORIZONTAL
  445.     5     X1=X+DX
  446.           Y1=Y-DY
  447.           X2=X-DX
  448.           Y2=Y+DY
  449.           GO TO 11
  450. C UPPER-LEFT
  451.     6     ASSIGN 12 TO IRET
  452.     7     X1=X+DX
  453.           Y1=Y-DY
  454.           X2=X
  455.           Y2=Y+DZ
  456.           GO TO 11
  457. C VERTICAL
  458.     8     X1=X
  459.           Y1=Y-DZ
  460.           X2=X
  461.           Y2=Y+DZ
  462.           GO TO 11
  463.     9     ASSIGN 12 TO IRET
  464. C UPPER-LEFT
  465.    10     X1=X-DX
  466.           Y1=Y+DY
  467.           X2=X
  468.           Y2=Y+DZ
  469. C TEST VISIBILITY OF THIS LINE SEGMENT.
  470.    11     IX=(X1-S(1))*RX
  471.           IY=MOD(IFIX((Y1-S(3))*RY-SLOPE*FLOAT(IX))+NY,NY)+1
  472.           IBIT=MOD(IX,NBPW)+1
  473.           IX=IX/NBPW+1
  474. C********* .AND. USED AS A MASKING OPERATOR ************** IAND KM Crennell
  475.           IV=IAND(IS2(IX,IY),MASK(IBIT))
  476. C IF EITHER END OF THE LINE IS AT A MARKED SPOT ON THE IMAGE
  477. C PLANE MODEL, THE LINE IS HIDDEN
  478.           IF(IV.NE.0) GO TO IRET,(6,9,12)
  479.           IX=(X2-S(1))*RX
  480.           IY=MOD(IFIX((Y2-S(3))*RY-SLOPE*FLOAT(IX))+NY,NY)+1
  481.           IBIT=MOD(IX,NBPW)+1
  482.           IX=IX/NBPW+1
  483. C******** .AND. USED AS A MASKING OPERATOR ******changed to IAND KM Crennell
  484.           IV=IAND(IS2(IX,IY),MASK(IBIT))
  485.           IF(IV.NE.0) GO TO IRET,(6,9,12)
  486. C*************** UNDEFINED EXTERNAL REFERENCE **************
  487. C SUBROUTINE LINE(X1,Y1,X2,Y2) IS ASSUMED TO DRAW A LINE
  488. C FROM (X1,Y1) TO (X2,Y2)
  489.           CALL CLINE(X1,Y1,X2,Y2)
  490.           GO TO IRET,(6,9,12)
  491.    12     CONTINUE
  492. C CODE THROUGH STATEMENT 13 CREATES AN APPROXIMATION OF
  493. C THE SILHOUETTE OF THE PART OF THE PICTURE JUST DRAWN BY
  494. C MARKING THE IMAGE PLANE MODEL WHERE THE OBJECT OCCURS.
  495.       DO 13 I=1,NV
  496.         DO 13 J=1,NW
  497.           IF(IOBJS(I,J).EQ.0) GO TO 13
  498.           IX=(ST1(I,J,1)-S(1))*RX+0.5
  499.           TWK=SLOPE*FLOAT(IX)-0.5
  500.           IY=MOD(IFIX((ST1(I,J,2)-S(3))*RY-TWK)+NY,NY)+1
  501.           IBIT=MOD(IX,NBPW)+1
  502.           IX=IX/NBPW+1
  503. C************ .OR. USED AS A MASKING OPERATOR *********made IOR KM Crennell
  504.           IS2(IX,IY)=IOR(IS2(IX,IY),MASK(IBIT))
  505.    13     CONTINUE
  506.       RETURN
  507.   100 IER=1
  508.       RETURN
  509.       END
  510.