C     STRCON
C
C     THIS PROGRAM CONTOURS FOURIER MAPS PRODUCING
C     AN UNDISTORTED MAP VIA THE VERSATEC MODEL 1100.
C     THE MAP CONTOURING PORTION IS BASED ON THE ORIGINAL
C      PROGRAM TRICON WRITTEN BY JOHN ROSENBERG.
C
C	MODIFIED EXTENSIVELY BY G.J.QUIGLEY, APRIL, 1978
C	FOR USE WITH VERSATEC PLOT 07 SOFTWARE
C       AND TO INCLUDE THE MODEL, VAN DER WALLS SURFACE AND GRID OPTIONS
C       AND FIX UP FOR NON ORTHOGANAL CELLS.
C       UPDATED NOV. 1979 G.J.Q.
C
C --------------------------------------------------------------
C
C	MODIFIED BY R.S. ZHOU, BOSTON COLLEGE, DEC, 1988
C	TO LINK TO GRAFIT'S PLOTTING LIBRARIES.
C
C***********************************************************************
C
	DIMENSION ID(3), A(3), ANG(3)
	DIMENSION XMIN2(3)
	INTEGER*4 PWIDE,CHAR,VTYP
	LOGICAL*1 ATMFIL
	COMMON/MAPDAT/  XMAP(150), YMAP(150), XMAPY(150)
	COMMON/STRO/ISTRO,STROX,STROY,SSTROX,SSTROY,MSEC,MSEC2,IINZ,NZMX
	COMMON /DAT/ TI(10), XMIN(3), XMAX(3), DX(3), XCHAR(3), ACHAR(3),
	1  IDIM(2),JDIM(2),NEXT,XSHX,XSHY,XSHZ,YSHY,YSHZ,ZSHZ,XMPMIN,
	2  XMPMAX, YMPMIN, YMPMAX, ISCALE, NSDONE, NIP, ZMAP, A,
	3  RMIN,RMAX,RINC,NSEC,ANG,RC,RMINN,RMAXN,RINCN
	COMMON /IOSET/ IO1,IO2,IO3,IO4,IO5,IO6,IO7,IO8
	COMMON /MOD/ AM(3,3,36),AT(3,36),XMIN3(3),XMAX3(3),
	1 DZ3,DZ3N,DZ2,DZ2N,DMX(4),CHAR(4,10),NSYM,PWIDE(4),ICHR(4),
	2 HT(4),HT2(4),SCALE,RAD(4,10),RDEL(4),IMASK(36,4),
	3 NSET(4,20), ATMFIL(36,4),VRAD(10),VTYP(10),IFVDW

     	INTEGER   CRT
     	LOGICAL   YES
     	CHARACTER MAPFIL*32, PARFIL*32, DESFIL*32
     	CHARACTER COMM*16, PLOTID*32
     	DIMENSION HPAR(20)
C
C  PARSE THE ACTIVATING DCL COMMAND LINE...
C
C     	CRT = IF_DCL_PUT( 'STRCON' )
C     	IF( CRT.LT.0 ) GOTO 888
     	CALL IFINIT_SESSION
     	IF( .NOT.IF_DCL_HARDCOPY() ) THEN
     	    CRT = IFINIT_ITM(0)
     	    CALL IFANSI_CLEAR(0)
     	ENDIF
     	XORG = 1.5
     	YORG = 1.5
     	PLOTID = ' '
C
C  GET THE MAP FILES...
C
     	CALL IF_DCL_FILE( 1, MAPFIL )
 4	IF( MAPFIL.NE.' ' ) THEN
     	    CALL IF_DCL_SETUP( PARFIL )
     	    IF( PARFIL.EQ.' ' ) THEN
     		PARFIL = MAPFIL
     		CALL IFMISC_NAMEXT( PARFIL, '.CTL' )
     		INQUIRE( FILE=PARFIL, EXIST=YES )
     		IF(.NOT.YES) STOP '%STRCON: missing setup file'
     	    ENDIF
     	    CALL IF_DCL_LAYOUT( DESFIL )
     	    IF( DESFIL.EQ.' ' ) THEN
     		DESFIL = MAPFIL
     		CALL IFMISC_NAMEXT( DESFIL, '.DES' )
     	    ENDIF
     	ELSE
 5     	    WRITE(6,5000)
5000	    FORMAT(' Density data file [.MAP] ',$)
     	    READ(5,5001,END=999) MAPFIL
5001	    FORMAT(A)
     	    IF( MAPFIL.EQ.' ' ) GOTO 5
     	    CALL IFMISC_DEFEXT( MAPFIL, '.MAP' )
     	    INQUIRE( FILE=MAPFIL, EXIST=YES )
     	    IF( .NOT.YES ) THEN
     		WRITE(6,6)
 6		FORMAT(' oops: Data file does not exist.')
     		GOTO 5
     	    ENDIF
     	    PARFIL = MAPFIL
     	    CALL IFMISC_NAMEXT( PARFIL, '.CTL' )
     	    INQUIRE( FILE=PARFIL, EXIST=YES )
     	    IF( YES ) GOTO 8
 7	    WRITE(6,5005)
5005	    FORMAT(' Setup file [.CTL] ',$)
     	    READ(5,5001,END=999) PARFIL
     	    IF( PARFIL.EQ.' ' ) GOTO 7
     	    CALL IFMISC_DEFEXT( PARFIL, '.CTL' )
     	    INQUIRE( FILE=PARFIL, EXIST=YES )
     	    IF( .NOT.YES ) THEN
     		WRITE(6,9)
 9		FORMAT(' oops: Setup file does not exist.')
     		GOTO 7
     	    ENDIF
 8	    DESFIL = MAPFIL
     	    CALL IFMISC_NAMEXT( DESFIL, '.DES' )
     	    XORG = 1.5
     	    YORG = 1.5
5101	    WRITE(6,5100)
5100	    FORMAT(' Plot offset (x,y): ',$)
     	    READ(5,5102,ERR=5101,END=999) ORGX,ORGY
5102	    FORMAT(2F10.0)
     	    IF( ORGX.GT.0.0 ) XORG = XORG + ORGX
     	    IF( ORGY.GT.0.0 ) YORG = YORG + ORGY
     	ENDIF
C
C  GET CONTOUR SETUP...
C
     	CALL STRCON_TRSTUP( PARFIL,
     +  XMIN2,MDIR,NDZ,RSCX,RSCY,XFRAM,YFRAM,IFNFLD,XFNFLD,IHEX)
C
C      leave 8 mms at the top for the label
C      automatically put the origin in the upper left corner
C      at last, begin contouring
C
	ZMAPS = XMIN(3)
	SC=1.0/SCALE
	SCX=1.0/(SCALE*RSCX)
	SCY=1.0/(SCALE*RSCY)
	SSTROX=DX(3)*STROX
	SSTROY=DX(3)*STROY
	STML0=YMPMIN
	STMT0=XMPMIN-8.0/SCALE
	STMT0=XMPMIN
	IF(ISTRO.NE.0)STML0=STML0-FLOAT(MSEC-1)*0.5*SSTROY
	IF(ISTRO.NE.0)STMT0=STMT0-FLOAT(MSEC-1)*0.5*SSTROX
	ASTROX=STROX
	ASTROY=STROY
C
C  START THE FRAME LOOP...
C
     	III = 0
 1234	III = III + 1
     	IF( III.LT.1 .OR. III.GT.NSEC ) THEN
     	    IF( IF_DCL_HARDCOPY() ) GOTO 4321
     	    CALL IFANSI_VIDEO(70)
     	    WRITE(6,1235)
 1235	    FORMAT(' oops: Frame number out of bound.')
     	    IF( III.LT.1 ) III = 1
     	    IF( III.GT.NSEC ) III = NSEC
     	    GOTO 444
     	ENDIF
     	IF( IF_DCL_HARDCOPY() ) WRITE(6,1236) III
 1236	FORMAT(' Contouring frame #',I2,'...',$)
	ZMAP = ZMAPS + MSEC2*DX(3)*(III-1)
	ZMAPO=ZMAP*ZSHZ
	ZMAPF=ZMAP/A(3)
	WRITE(IO6, 230) XCHAR(3), ZMAPO, ZMAPF
 230	FORMAT('0BEGINNING SECTION ', A1, ' = ', F7.3,'A (',F6.3,')')
C
C  LABEL SECTION AND DRAW GRID...
C
	STROSN=0
	DO 260 JSTRO=1,2
	IF(ISTRO.EQ.0.AND.JSTRO.EQ.2)GO TO 260
C
C  INITIALIZE THE PLOTTER...
C
     	IF( JSTRO.EQ.1 ) THEN
     	    CALL IFINIT_LINE(3)
     	    CALL IFINIT_PLOT( XORG, YORG, MAPFIL )
     	    IF( IFCTRL_C() ) GOTO 333
C
C  OPEN A GRAPHICS SEGMENT TO SAVE THE PLOT...
C
     	    IF( CRT.EQ.2 ) THEN
     		CALL IFSEG_PIVOT( 0.0, 0.0 )
     		CALL IFSEG_OPEN( 15 )
     	    ENDIF
     	    CALL IFINIT_DIF(0)
     	    CALL IFSPEN(7)
     	    COMM = 'CURRENT_FRAME:'
     	    ENCODE( 2, 240, COMM(15:16) ) III
  240	    FORMAT(I2)
     	    CALL IFTEXT( 0.1-XORG, 7.0-YORG, COMM, 0.1, 0.1, 0.0 )
     	    CALL IFINIT_DIF(1)
     	    CALL IFSPEN(1)
     	ENDIF
	IF(ISTRO.NE.0)STROSN=(-1.0)**JSTRO
	STROX=ASTROX*STROSN
	STROY=ASTROY*STROSN
     	IF( IFCTRL_C() ) GOTO 333
     	CALL IFSPEN(1)
	CALL STRCON_FRAME( SC,RSCX,RSCY,XFRAM,YFRAM,IFNFLD,XFNFLD)
     	IF( IFCTRL_C() ) GOTO 333
	IF(RMAX.NE.0.0.OR.RMAXN.NE.0.0) THEN
	   OPEN(UNIT=IO1,FILE=MAPFIL,TYPE='OLD',FORM='UNFORMATTED')
	   INZ=(IINZ+MSEC2*(III-1))*NDZ
	   IF(INZ.GT.0)THEN
	      DO II=1,INZ
	         READ(IO1)
	      END DO
	   END IF
	END IF
	CALL STRCON_OFFSET(STMT0,SCX,STML0,SCY)
     	IF( IFCTRL_C() ) GOTO 333
	CALL STRCON_FACE(ISTRO,JSTRO)
     	IF( IFCTRL_C() ) GOTO 333
C
C  START THE SECTION STACKING LOOP...
C
	DO 250 I=1,MSEC
     	IF( IFCTRL_C() ) GOTO 333
	STML=STML0-0.5*FLOAT(2*I-1-MSEC)*SSTROY*STROSN
	STMT=STMT0-0.5*FLOAT(2*I-1-MSEC)*SSTROX*STROSN
	CALL STRCON_OFFSET(STMT,SCX,STML,SCY)
	IF((IINZ+MSEC2*(III-1)+I).LE.0.OR.
	1  (IINZ+MSEC2*(III-1)+I).GT.NZMX)GO TO 250
	IF(RMAX.EQ.0.0.AND.RMAXN.EQ.0.0)GO TO 250
C
C  DEFINE THE FLOATING GRID...
C
	CALL STRCON_GRID( XMAP, YMAP, XMAPY )
     	IF( IFCTRL_C() ) GOTO 333
C
C  GET THE SECTION...
C
	IF(MDIR.GT.0) CALL STRCON_SECTN (XMIN2)
     	IF( IFCTRL_C() ) GOTO 333
	IF(MDIR.LT.0) CALL STRCON_SECTNN(XMIN2)
     	IF( IFCTRL_C() ) GOTO 333
C
C  DO THE ACTUAL CONTOURING...
C
     	IF(RMAX.NE.0.0) THEN
     	CALL IFSPEN(3)
     	IF(IHEX.EQ.0) CALL STRCON_TACON(JDIM(1),JDIM(2),RMIN,RMAX,RINC)
     	IF(IHEX.NE.0) CALL STRCON_HACON(JDIM(1),JDIM(2),RMIN,RMAX,RINC)
     	ENDIF
     	IF( IFCTRL_C() ) GOTO 333
	IF(NONEG.NE.1.AND.RMAXN.NE.0.0) THEN
     	CALL IFSPEN(4)
	IF(IHEX.EQ.0) CALL STRCON_TACOF(JDIM(1),JDIM(2),RMINN,RMAXN,
     +  RINCN)
	IF(IHEX.NE.0) CALL STRCON_HACOF(JDIM(1),JDIM(2),RMINN,RMAXN,
     +  RINCN)
     	ENDIF
     	CLOSE(UNIT=IO3)
     	IF( IFCTRL_C() ) GOTO 333
 250	CONTINUE
	IF(RMAX.NE.0.0.OR.RMAXN.NE.0.0)CLOSE(UNIT=IO1)
	CALL STRCON_OFFSET(STMT0,SCX,STML0,SCY)
     	IF( IFCTRL_C() ) GOTO 333
     	CALL IFSPEN(2)
 	CALL STRCON_MODEL(STROX,STROY)
     	IF( IFCTRL_C() ) GOTO 333
	IF(IFVDW.NE.0) THEN
     	CALL IFSPEN(5)
     	CALL STRCON_VDWDRW(SSTROX,SSTROY,STML0,STMT0,MSEC,STROSN,SCX,SCY)
     	IF( IFCTRL_C() ) GOTO 333
	END IF
	NSDONE = NSDONE + 1
	IF(IFNFLD.NE.0) THEN
	   XXXX = XFNFLD * RSCX
	   CALL STRCON_QPLOT( XXXX, 0.0, -3)
	END IF
 260	CONTINUE
	WRITE(IO6, 231)
 231	FORMAT('0SECTION COMPLETED.')
C
C  CLOSE THE FRAME...
C
     	IF( CRT.EQ.2 ) CALL IFSEG_OPEN( 0 ) 	!CLOSE THE SEGMENT
    	IF( DESFIL.NE.' ' ) CALL LAYOUT_AUX(0,XORG,YORG,DESFIL)
  333	CALL IFINIT_PLOT_END
     	IF( IF_DCL_HARDCOPY() ) GOTO 666
     	IF( CRT.EQ.1 ) THEN
     	    CALL IFANSI_CURSOR(1,1)
     	    CALL IFANSI_CLEAR(6)
     	ENDIF
C
C  PAUSE FOR INTERACTIVE INPUT....
C
  444	CALL IFANSI_VIDEO(47)
     	WRITE(6,1000)
 1000	FORMAT(' STRCON> Please press a desired function key:')
  445	WRITE(6,1001)
 1001	FORMAT(' STRCON> ?>Help 0>Exit 1>Hardcopy 2>Layout',
     +' 3>Newfile 4>Down_ 5>Up_frame ...')
     	CALL IFANSI_CURSOR(0,8)
     	CALL IFQIOW_KEY( KEY )

     	IF( KEY.EQ.26 .OR. KEY.EQ.48 .OR. KEY.EQ.69 ) GOTO 4321
     	IF( KEY.EQ.49 .OR. KEY.EQ.72 ) THEN		!hardcopy
     	    CALL IFCOPY( HPAR(1), PLOTID )
     	    CALL IFINIT_LINE( -3 )
     	ELSE IF( KEY.EQ.-80 ) THEN			!hardocpy
     	    CALL IFANSI_CURSOR(0,14)
     	    CALL IFANSI_LEDS(1)
     	    CALL IFFEED( HPAR(1), PLOTID, ' ' )
     	    CALL IFANSI_LEDS(0)
     	ELSE IF( KEY.EQ.50 .OR. KEY.EQ.76 ) THEN	!layout
     	    IF( DESFIL.EQ.' ' ) THEN
     		DESFIL = MAPFIL
     		CALL IFMISC_NAMEXT( DESFIL, '.DES' )
     	    ENDIF
     	    CALL LAYOUT_AUX( 2, XORG, YORG, DESFIL )
     	ELSE IF( KEY.EQ.51 .OR. KEY.EQ.-82 ) THEN	!newfile
     	    MAPFIL = ' '
    	    GOTO 4
     	ELSE IF( KEY.EQ.52 .OR. KEY.EQ.68 ) THEN	!down_frame
    	    GOTO 1234
     	ELSE IF( KEY.EQ.53 .OR. KEY.EQ.85 ) THEN	!up_frame
     	    III = III-2
     	    GOTO 1234
     	ELSE IF( KEY.EQ.54 ) THEN	!redraw the same frame
     	    III = III-1
     	    GOTO 1234
     	ELSE IF( KEY.EQ.13 .OR. KEY.EQ.63
     +  .OR. KEY.EQ.32 .OR. KEY.EQ.-81 ) THEN	!help
     	    CALL IFMISC_HLP( 'STRCON', '$ MENU', -2, 0 )
     	ELSE IF( KEY.EQ.36 .OR. KEY.EQ.-83 ) THEN	!dcl prompt
     	    CALL IFMISC_SPAWN( '$', 0 )
     	    CALL IFINIT_LINE( -3 )
     	ELSE IF( KEY.EQ.64 ) THEN
554	    WRITE(6,555)
555	    FORMAT(' Which frame: ',$)
     	    READ(5,556,ERR=554,END=444) KEY
556	    FORMAT(I)
     	    IF( KEY.LT.1 .OR. KEY.GT.NSEC ) THEN
     		CALL IFANSI_VIDEO(70)
     		WRITE(6,557)
557		FORMAT('_oops: non-existing frame.')
     		GOTO 554
     	    ENDIF
     	    III = KEY - 1
     	    GOTO 1234
     	ELSE
     	    CALL IFMISC_OOPS( 5, 0, 0, ' ' )
    	ENDIF
     	GOTO 445
C
C  HARDCOPY THE CURRENT FRAME...
C
 666	CALL IFCOPY( HPAR(1), PLOTID )
     	GOTO 1234			!GET THE NEXT FRAME
4321	CLOSE(UNIT=IO6)
999	CALL IFINIT_EXIT
888	CALL EXIT
	END

	SUBROUTINE STRCON_QUIT(NSDONE)
	COMMON /IOSET/IO1,IO2,IO3,IO4,IO5,IO6,IO7,IO8
	WRITE(IO6, 100) NSDONE
 100	FORMAT('0HALTING AFTER COMPLETING ', I2, ' SECTIONS.')
	IF(NSDONE .EQ. 0) STOP 'ABNORMAL TERMINATION TO STRCON.'
	CLOSE(UNIT=IO6)
	STOP 'STRCON NORMAL TERMINATION.'
	END

	SUBROUTINE STRCON_ROT(X,Y,Z,ISORT)
C
C
	IF(ISORT.EQ.0) RETURN
	GO TO (1,4,2,5,3,6) ISORT
	STOP 'NO SUCH SORT OPTION'
C -XYZ-
1	RETURN
C -XZY-
2	T=Y
	Y=Z
	Z=T
	RETURN
C -YZX-
3	T=X
	X=Y
	Y=Z
	Z=T
	RETURN
C -YXZ-
4	T=X
	X=Y
	Y=T
	RETURN
C -ZXY -
5	T=X
	X=Z
	Z=Y
	Y=T
	RETURN
6	T=X
C -ZYX-
	X=Z
	Z=T
	RETURN
	END

	SUBROUTINE STRCON_GRID( XMAP, YMAP, XMAPY )
C
C      DEFINES A FLOATING GRID COORECTLY POSITIONED FOR THIS SECTION
C
	COMMON /DAT/ TI(10), XMIN(3), XMAX(3), DX(3), XCHAR(3), ACHAR(3),
	1  IDIM(2),JDIM(2),NEXT,XSHX,XSHY,XSHZ,YSHY,YSHZ,ZSHZ,XMPMIN,
	2  XMPMAX, YMPMIN, YMPMAX, ISCALE, NSDONE, NIP, ZMAP, A(3),
	3  RMIN, RMAX, RINC, NSEC, ANG(3), RC
	DIMENSION XMAP(150), YMAP(150), XMAPY(150)
	XORIG = XMIN(1) + XSHY*XMIN(2) + XSHZ*ZMAP
	YORIG = YSHY*XMIN(2) + YSHZ*ZMAP
	JDL1 = (XORIG - XMPMIN)/(DX(1)*XSHX) + 0.50
	JDH1 = (XMPMAX - XORIG)/(DX(1)*XSHX) + 0.50
	JDL2 = (YORIG - YMPMIN)/(YSHY*DX(2)) + 0.50
	JDH2 = (YMPMAX - YORIG)/(YSHY*DX(2)) + 0.50
	JDIM1 = JDL1 + JDH1 + 1
	JDIM2 = JDL2 + JDH2 + 1
	XYOR=YMPMIN
	IF(XSHY.LT.0.0)XYOR=YMPMAX
	DO 1 IND = 1, JDIM1
	X = IND - JDL1 - 1
	XMAP(IND) = XORIG + X*DX(1)*XSHX
 1	CONTINUE
	DO 2 IND = 1, JDIM2
	Y = IND - JDL2 - 1
	YMAP(IND) = YORIG + Y*(YSHY*DX(2))
	Y1=(YMAP(IND)-XYOR)*XSHY/(YSHY*DX(1))
	IY = IFIX(10000.0 + Y1) - 10000
	Y1=Y1-IY
	IF(Y1.GT.0.51)Y1=Y1-1.0
	XMAPY(IND)=Y1*DX(1)
 2	CONTINUE
	RETURN
	END

	SUBROUTINE STRCON_SECTN(XMIN2)
C
C      GETS THE SECTION, SHIFTS COORDINATES, AND EXTENDS THE MAP
C
	COMMON/MAPDAT/  XMAP(150), YMAP(150), XMAPY(150)
	COMMON /DAT/ TI(10), XMIN(3), XMAX(3), DX(3), XCHAR(3), ACHAR(3),
	1 IDIM1,IDIM2,JDIM1,JDIM2,NEXT,XSHX,XSHY,XSHZ,YSHY,YSHZ,ZSHZ,
	2  XMPMIN,XMPMAX,YMPMIN,YMPMAX,ISCALE,NSDONE,NIP,ZMAP,A(3),
	3  RMIN, RMAX, RINC, NSEC, ANG(3), RC
	COMMON /IOSET/ IO1,IO2,IO3,IO4,IO5,IO6,IO7,IO8
	LOGICAL INSIDE
	INTEGER*2 MAP(100),RMAP(109,109),RB
	DIMENSION XMIN2(3)
	DIMENSION COORD(2, 150), XX(2)
	OPEN(UNIT=IO3,STATUS='SCRATCH',FORM='UNFORMATTED')
C      COORD(1, I) IS THE X COORDINATE OF GRID POINT I, J.
C      COORD(2, J) IS THE Y COORDINATE.
C      GENERATE UNSHIFTED COORDINATES
	DO 1 I = 1, 150
	XX(1) = I - 1
	COORD(1, I) = XMIN2(1) + DX(1)*XX(1) + XSHY*XMIN2(2)
	COORD(2, I) = YSHY*XMIN2(2) + DX(2)*XX(1)*YSHY
 1	CONTINUE
C      CLEAR CONTOURING ARRAY
	RB= 32767
	DO 4 I = 1, JDIM1
	DO 4 J = 1, JDIM2
	RMAP( I,J) = RB
 4	CONTINUE
C      GENERATE CONTOURING ARRAY
	DO 16 I = 1, IDIM1
	READ(IO1,END=1000,ERR=2000)(MAP(I1),I1=1,IDIM2)
	DO 6 J = 1, IDIM2
	DO 6 K = 1, 3
	DO 6 L = 1, 3
	GO TO (7, 8, 9), K
 7	CONTINUE
	IF((NEXT .EQ. 0) .OR. (NEXT .EQ. 2)) GO TO 6
	XSHIF = -A(1)
	GO TO 10
 8	CONTINUE
	IF((NEXT .EQ. 0) .OR. (NEXT .EQ. 2)) GO TO 6
	XSHIF = A(1)
	GO TO 10
 9	CONTINUE
	XSHIF = 0
 10	CONTINUE
	XSHIF = XSHIF + COORD(1, I)
	GO TO (11, 12, 13) , L
 11	CONTINUE
	IF((NEXT .EQ. 0) .OR. (NEXT .EQ. 1)) GO TO 6
	YSHIF = -A(2)
	GO TO 14
 12	CONTINUE
	IF((NEXT .EQ. 0) .OR. (NEXT .EQ. 1)) GO TO 6
	YSHIF = A(2)
	GO TO 14
 13	CONTINUE
	YSHIF = 0
 14	CONTINUE
	YSHIF = YSHIF*YSHY + COORD(2, J)
	XSHIF = XSHIF + XSHY*YSHIF/YSHY + XSHZ*ZMAP
	YSHIF = YSHIF + YSHZ*ZMAP
	CALL STRCON_INDEX( XSHIF, YSHIF, IX, IY, INSIDE)
	IF(.NOT. INSIDE) GO TO 6
	RMAP(IX,IY) = MAP(J)
 6	CONTINUE
 16	CONTINUE
	DO 26 J=1,JDIM2
	WRITE(IO3)(RMAP(I,J),I=1,JDIM1)
 26	CONTINUE
	REWIND IO3
	RETURN
 1000	WRITE(IO6, 1001)
 1001	FORMAT('0REACHED THE END OF THE MAP DATASET.')
	CALL STRCON_QUIT(NSDONE)
	RETURN
 2000	NN = NSDONE + 1
	WRITE(IO6, 2001) NN
 2001	FORMAT('0ENCOUNTERED AN ERROR WHILE TRYING TO READ SECTION ', I2)
	CALL STRCON_QUIT(NSDONE)
	RETURN
	END

	SUBROUTINE STRCON_SECTNN(XMIN2)
C
C      GETS THE SECTION, SHIFTS COORDINATES, AND EXTENDS THE MAP
C
	COMMON/MAPDAT/  XMAP(150), YMAP(150), XMAPY(150)
	COMMON /DAT/ TI(10), XMIN(3), XMAX(3), DX(3), XCHAR(3), ACHAR(3),
	1  IDIM1,IDIM2,JDIM1,JDIM2,NEXT,XSHX,XSHY,XSHZ,YSHY,YSHZ,ZSHZ,
	2 XMPMIN,XMPMAX,YMPMIN,YMPMAX,ISCALE,NSDONE,NIP,ZMAP,A(3),
	3  RMIN, RMAX, RINC, NSEC, ANG(3), RC
	COMMON /IOSET/ IO1,IO2,IO3,IO4,IO5,IO6,IO7,IO8
	LOGICAL INSIDE
	INTEGER*2 MAP(100),RMAP(109,109),RB
	DIMENSION XMIN2(3)
	DIMENSION COORD(2, 150), XX(2)
	OPEN(UNIT=IO3,STATUS='SCRATCH',FORM='UNFORMATTED')
C      COORD(1, I) IS THE X COORDINATE OF GRID POINT I, J.
C      COORD(2, J) IS THE Y COORDINATE.
C      GENERATE UNSHIFTED COORDINATES
	DO 1 I = 1, 150
	XX(1) = I - 1
	COORD(1, I) = XMIN2(1) + DX(1)*XX(1) + XSHY*XMIN2(2)
	COORD(2, I) = YSHY*XMIN2(2) + DX(2)*XX(1)*YSHY
 1	CONTINUE
C      CLEAR CONTOURING ARRAY
	RB= 32767
	DO 4 I = 1, JDIM1
	DO 4 J = 1, JDIM2
	RMAP( I,J) = RB
 4	CONTINUE
C      GENERATE CONTOURING ARRAY
	DO 16 J = 1, IDIM2
	READ(IO1,END=1000,ERR=2000)(MAP(I1),I1=1,IDIM1)
	DO 6 I = 1, IDIM1
	DO 6 K = 1, 3
	DO 6 L = 1, 3
	GO TO (7, 8, 9), K
 7	CONTINUE
	IF((NEXT .EQ. 0) .OR. (NEXT .EQ. 2)) GO TO 6
	XSHIF = -A(1)
	GO TO 10
 8	CONTINUE
	IF((NEXT .EQ. 0) .OR. (NEXT .EQ. 2)) GO TO 6
	XSHIF = A(1)
	GO TO 10
 9	CONTINUE
	XSHIF = 0
 10	CONTINUE
	XSHIF = XSHIF + COORD(1, I)
	GO TO (11, 12, 13) , L
 11	CONTINUE
	IF((NEXT .EQ. 0) .OR. (NEXT .EQ. 1)) GO TO 6
	YSHIF = -A(2)
	GO TO 14
 12	CONTINUE
	IF((NEXT .EQ. 0) .OR. (NEXT .EQ. 1)) GO TO 6
	YSHIF = A(2)
	GO TO 14
 13	CONTINUE
	YSHIF = 0
 14	CONTINUE
	YSHIF = YSHIF*YSHY + COORD(2, J)
	XSHIF = XSHIF + XSHY*YSHIF/YSHY + XSHZ*ZMAP
	YSHIF = YSHIF + YSHZ*ZMAP
	CALL STRCON_INDEX( XSHIF, YSHIF, IX, IY, INSIDE)
	IF(.NOT. INSIDE) GO TO 6
	RMAP(IX,IY) = MAP(I)
 6	CONTINUE
 16	CONTINUE
	DO 26 J=1,JDIM2
	WRITE(IO3)(RMAP(I,J),I=1,JDIM1)
 26	CONTINUE
	REWIND IO3
	RETURN
 1000	WRITE(IO6, 1001)
 1001	FORMAT('0REACHED THE END OF THE MAP DATASET.')
	CALL STRCON_QUIT(NSDONE)
	RETURN
 2000	NN = NSDONE + 1
	WRITE(IO6, 2001) NN
 2001	FORMAT('0ENCOUNTERED AN ERROR WHILE TRYING TO READ SECTION ', I2)
	CALL STRCON_QUIT(NSDONE)
	RETURN
	END

	SUBROUTINE STRCON_FRAME(SC,RSCX,RSCY,XFRAM,YFRAM,IFNFLD,XFNFLD)

	CALL STRCON_QPLOT(0.0,0.0,+3)
	XX=XFRAM*RSCX
	XXX=XFNFLD*RSCX
	YY=YFRAM*RSCY
	CALL STRCON_QPLOT(0.0,YY,+2)
	CALL STRCON_QPLOT(XX,YY,+2)
	CALL STRCON_QPLOT(XX,0.0,+2)
	CALL STRCON_QPLOT(0.0,0.0,+2)
	IF(IFNFLD.NE.0.0) CALL STRCON_QPLOT(XXX,0.0,+2)
	RETURN
	END

	SUBROUTINE STRCON_FACE(ISTRO,JSTRO)
C
C      Labels the section, etc.
C
	COMMON /DAT/ TI(10), XMIN(3), XMAX(3), DX(3), XCHAR(3), ACHAR(3),
	1  IDIM(2),JDIM(2),NEXT,XSHX,XSHY,XSHZ,YSHY,YSHZ,ZSHZ,XMPMIN,
	2  XMPMAX, YMPMIN, YMPMAX, ISCALE, NSDONE, NIP, ZMAP, A(3),
	3  RMIN, RMAX, RINC, NSEC, ANG(3), RC
	COMMON /IOSET/ IO1,IO2,IO3,IO4,IO5,IO6,IO7,IO8
	REAL*4 IX1,IY1,IX2,IY2,IX3,IY3,IX4,IY4
	REAL*4 IX5,IY5,IX6,IY6,IX7,IY7,IX8,IY8
     	CHARACTER TITLE*40, SECT*40
	LOGICAL*1 STRCHR(5),STROC
	DATA STRCHR/' ','T','B','L','R'/
     	EQUIVALENCE (TITLE(1:1),TI(1))
C
C      Title the section
C
     	IF( TITLE.NE.' ' ) THEN
     	    CALL IFSPEN(3)
     	    CALL IFTEXT( 0.0, -0.6, TITLE, 0.12, 0.12, 0.0 )
     	ENDIF
	K=1+JSTRO+2*(ISTRO-1)
	IF(ISTRO.EQ.0)K=1
	STROC=STRCHR(K)
	ISC=IABS(ISCALE)
     	ENCODE(40,1,SECT) ISC,XCHAR(3),ZMAP,ACHAR(2),ACHAR(1),STROC
  1	FORMAT('G=',I2,'A ',A1,'=',F6.2,'A, ',A1,
     +  '-HORI,',A1,'-VERT, S.P.-',A1,'.')
     	CALL IFSPEN(1)
	CALL IFTEXT( 0.0, -0.2, SECT, 0.1, 0.1, 0.0 )
C
C      outline the plotting area
C
	IX1 = XMPMIN
	IX2 = XMPMAX
	IY1 = YMPMIN
	IY2 = YMPMAX
	CALL STRCON_LINEV(IX1, IY1, IX1, IY2)
	CALL STRCON_LINEV(IX2, IY1, IX2, IY2)
	CALL STRCON_LINEV(IX1, IY1, IX2, IY1)
	CALL STRCON_LINEV(IX1, IY2, IX2, IY2)
C
C      place scaled grid
C
	IF(ISCALE.EQ.0)GO TO 6
	SC=IABS(ISCALE)
	SC1=SC/10.0
	NTKXI=IFIX(XMPMIN/SC + 1000.0)- 1000
	NTKYI=IFIX(YMPMIN/SC + 1000.0)- 1000
	NTKX=(XMPMAX-SC*NTKXI)/SC
	NTKY=(YMPMAX-SC*NTKYI)/SC
	IF(ISCALE.LT.0)GO TO 11
	NTKXI1=IFIX(XMPMIN/SC1 + 1000.0)- 1000
	NTKYI1=IFIX(YMPMIN/SC1 + 1000.0)- 1000
	NTKX1=(XMPMAX-SC1*NTKXI1)/SC1
	NTKY1=(YMPMAX-SC1*NTKYI1)/SC1
	DO 4 I = 1, NTKX
	   X = I + NTKXI
	   X = SC*X
	   DO 8 J=1,NTKY1
	      Y=J+NTKYI1
	      Y=SC1*Y
	      CALL STRCON_QPLOT(X,Y,+13)
	      CALL STRCON_QPLOT(X,Y,+12)
 8	   CONTINUE
 4	CONTINUE
	DO 5 I = 1, NTKY
	   Y = I + NTKYI
	   Y = SC*Y
	   DO 10 J=1,NTKX1
	      X=J+NTKXI1
	      X=X*SC1
	      CALL STRCON_QPLOT(X,Y,+13)
	      CALL STRCON_QPLOT(X,Y,+12)
 10	   CONTINUE
 5	CONTINUE
	GO TO 6
 11	CONTINUE
	DO 12 I=1,NTKX
	   X=I+NTKXI
	   X=X*SC
	   DO 12 J=1,NTKY
	      Y=SC*(J+NTKYI)
	      CALL STRCON_MARK(X,Y,0.08,9)
 12	   CONTINUE
 6	CONTINUE
C
C      outline unit cell section
C
	X = XMPMIN + XSHY*YMPMIN/YSHY + XSHZ*ZMAP
	IX5 = X
	X = XMPMAX + XSHY*YMPMIN/YSHY + XSHZ*ZMAP
	IX6 = X
	X = XMPMIN + XSHY*YMPMAX/YSHY + XSHZ*ZMAP
	IX7 = X
	X = XMPMAX + XSHY*YMPMAX/YSHY + XSHZ*ZMAP
	IX8 = X
	IF(IX5.GT.IX7)IX5=IX7
	IF(IX6.LT.IX8)IX6=IX8
	AA=A(1)/2.0
	N1=IFIX(IX5/AA+10.0)-10
	N2=IFIX(IX6/AA+10.0)-10
	N2=N2-N1+1
	N1=N1
C	write( 6, *) n1, n2, xmpmin, xmpmax, ympmin, ympmax
C	write( 6, *) xshy, yshy, xshz, zmap, aa
C	write( 6, *) ix5, ix6, ix7, ix8
	IF(N2.EQ.0)GO TO 22
	DO 22 I=1,N2
	   IY5=YMPMIN
	   IY6=YMPMAX
	   IX7=(I+N1)*AA+XSHY*YMPMIN/YSHY + XSHZ*ZMAP
	   IX8=(I+N1)*AA+XSHY*YMPMAX/YSHY + XSHZ*ZMAP
	   IF(IX7.LT.XMPMIN.AND.IX8.LT.XMPMIN)GO TO 22
	   IF(IX7.GT.XMPMAX.AND.IX8.GT.XMPMAX)GO TO 22
	   IF(IX7.LT.XMPMIN.OR. IX7.GT.XMPMAX)GO TO 24
	   IF(IX8.LT.XMPMIN.OR. IX8.GT.XMPMAX)GO TO 24
	   GO TO 20
 24	   IF(IX7.EQ.IX8)GO TO 22
 	   A1=(YMPMAX-YMPMIN)/(IX8-IX7)
	   B1=YMPMAX-A1*IX8
	   IF(IX7.LT.XMPMIN)IY5=A1*XMPMIN+B1
	   IF(IX7.LT.XMPMIN)IX7=XMPMIN
	   IF(IX7.GT.XMPMAX)IY5=A1*XMPMAX+B1
	   IF(IX7.GT.XMPMAX)IX7=XMPMAX
	   IF(IX8.LT.XMPMIN)IY6=A1*XMPMIN+B1
	   IF(IX8.LT.XMPMIN)IX8=XMPMIN
	   IF(IX8.GT.XMPMAX)IY6=A1*XMPMAX+B1
	   IF(IX8.GT.XMPMAX)IX8=XMPMAX
 20	   CONTINUE
C	   write( 6, *) xmpmin, xmpmax, a1, b1, aa
C	   write( 6, *) ix7, iy5, ix8, iy6
	   CALL STRCON_LINEV(IX7,IY5,IX8,IY6)
 22	CONTINUE
	Y = YMPMIN + YSHZ*ZMAP
	IY5 = Y
	Y = YMPMAX + YSHZ*ZMAP
	IY6 = Y
	IX5=XMPMIN
	IX6=XMPMAX
	AA = YSHY*A(2)/2.0
	N1=IFIX(IY5/AA+10.0)-10
	N2=IFIX(IY6/AA+10.0)-10
	N2=N2-N1
	N1=N1
C	write( 6, *) n1, n2, xmpmin, xmpmax, ympmin, ympmax
C	write( 6, *) xshy, yshy, xshz, zmap, aa
C	write( 6, *) ix5, ix6, ix7, ix8
	IF(N2.EQ.0)GO TO 32
	   DO 30 I=1,N2
	      Y=(I+N1)*AA+YSHZ*ZMAP
	      IY5=Y
C	      write( 6, *) ix5, iy5, ix6, iy5
	      CALL STRCON_LINEV(IX5, IY5, IX6, IY5)
 30	   CONTINUE
 32	CONTINUE
	RETURN
	END

	SUBROUTINE STRCON_INDEX( X, Y, IX, IY, INSIDE)
C
C      GETS INDICES FOR EXTENDED MAP ARRAY
C
	COMMON/MAPDAT/  XMAP(150), YMAP(150), XORY(150)
	COMMON /DAT/ TI(10),XMIN(3),XMAX(3),DX(3),XCHAR(3),ACHAR(3),
	1  IDIM(2),JDIM(2),NEXT,XSHX,XSHY,XSHZ,YSHY,YSHZ,ZSHZ,XMPMIN,
	2  XMPMAX, YMPMIN, YMPMAX, ISCALE, NSDONE, NIP, ZMAP, A(3),
	3  RMIN, RMAX, RINC, NSEC, ANG(3), RC
	LOGICAL INSIDE
	INSIDE = .FALSE.
	XORIGN = XMAP(1)
	YORIGN = YMAP(1)
	S = (Y - YORIGN)/(DX(2)*YSHY) + 1.
	IY = S
	YI = IY
	IF((S - YI) .GE. 0.5) IY = IY + 1
	IF((IY.LE.0).OR.(IY.GT.JDIM(2)))RETURN
	R = (X - XORIGN - XORY(IY))/DX(1) + 1.0
	IX = R
	XI = IX
	IF((R - XI) .GE. 0.51) IX = IX + 1
	IF((IX .GT. 0) .AND. (IX .LE. JDIM(1)) .AND. (IY .GT. 0)
	1  .AND. (IY .LE. JDIM(2))) INSIDE = .TRUE.
	RETURN
	END

	SUBROUTINE STRCON_LINEV(X1,Y1,X2,Y2)
C
C	SUBROUTINE TO DRAW A LINE FROM X1,Y1 TO X2,Y2.
C
	CALL STRCON_QPLOT(X1,Y1,+13)
	CALL STRCON_QPLOT(X2,Y2,+12)
	RETURN
	END

        SUBROUTINE STRCON_IROT(H,K,L,ISORT)
	INTEGER H,K,L
C
	IF(ISORT.EQ.0) RETURN
	GO TO (1,4,2,5,3,6) ISORT
	STOP 'NO SUCH SORT OPTION'
1	RETURN
2	I=K
	K=L
	L=I
	RETURN
 5	I=H
	H=K
	K=L
	L=I
	RETURN
4	I=H
	H=K
	K=I
	RETURN
3	I=H
	H=L
	L=K
	K=I
	RETURN
6	I=H
	H=L
	L=I
	RETURN
	END

	SUBROUTINE STRCON_TRSTUP( PARFIL,
     +  XMIN2,MDIR,NDZ,RSCX,RSCY,XFRAM,YFRAM,IFNFLD,XFNFLD,IHEX)
C
C	THIS ROUTINE GETS CONTOUR PARAMETERS.
C
	COMMON /DAT/ TI(10), XMIN(3), XMAX(3), DX(3), XCHAR(3), ACHAR(3),
	1  IDIM(2),JDIM(2),NEXT,XSHX,XSHY,XSHZ,YSHY,YSHZ,ZSHZ,XMPMIN,
	2  XMPMAX, YMPMIN, YMPMAX, ISCALE, NSDONE, NIP, ZMAP, A,
	3  RMIN,RMAX,RINC,NSEC,ANG,RC,RMINN,RMAXN,RINCN
	COMMON/STRO/ISTRO,STROX,STROY,SSTROX,SSTROY,MSEC,MSEC2,IINZ,NZMX
	COMMON /IOSET/ IO1,IO2,IO3,IO4,IO5,IO6,IO7,IO8
	COMMON /MOD/ AM(3,3,36),AT(3,36),XMIN3(3),XMAX3(3),
	1 DZ3,DZ3N,DZ2,DZ2N,DMX(4),CHR(4,10),NSYM,PWIDE(4),ICHR(4),
	2 HT(4),HT2(4),SCALE,RAD(4,10),RDEL(4),IMASK(36,4),
	3 NSET(4,20),ATMFIL(36,4),VRAD(10),VTYP(10),IFVDW
	DIMENSION ID(3), A(3), ANG(3)
	DIMENSION XMIN2(3), XMAX2(3)
	DIMENSION XMIN1(3),XMAX1(3)
	REAL*4 ML,MT
     	CHARACTER LOGFIL*32, PARFIL*(*)
	INTEGER*4 PWIDE,CHR,VTYP
	DIMENSION CHAR(6)
	DATA CHAR/'X', 'Y', 'Z', 'A', 'B', 'C'/
	DGTRD=0.017453293
	DO 30 I = 1, 3
	J = I + 3
	XCHAR(I) = CHAR(I)
 30	ACHAR(I) = CHAR(J)
	NIP = 12
	IO1 = 4
	IO2 = 2
	IO3 = 7
	IO4 = 4
	IO5 = 12
	IO6 = 13
     	LOGFIL = PARFIL
     	CALL IFMISC_NAMEXT( LOGFIL, '.LOG' )
     	IF( IF_DCL_HARDCOPY() ) WRITE(6,1) LOGFIL
 1	FORMAT(/' *** Run log file: ',A32,/)
	OPEN(UNIT=IO6,FILE=LOGFIL,FORM='FORMATTED',STATUS='UNKNOWN')
	RC = 3.141592653/180.0
	OPEN(UNIT=IO2,FILE=PARFIL,TYPE='OLD',FORM='FORMATTED')
	READ (IO2, 100) TI
 100	FORMAT(10A4)
	WRITE(IO6, 101) TI
 101	FORMAT(' ', 10A4)
C
C      get e.d. map parameters
C
	READ(IO2, 102) A, ANG
 102	FORMAT(3F10.4, 3F10.6)
	WRITE(IO6, 103) A, ANG
 103	FORMAT('0A = ', F10.4, 10X, 'B = ', F10.4, 10X, 'C = ', F10.4/
	1  '0ALPHA = ', F10.6, 10X, 'BETA = ', F10.6, 10X, 'GAMMA = ',
	2  F10.6)
	READ(IO2,1103)(XMIN1(I),XMAX1(I),I=1,3),MDIR
 1103	FORMAT(6F10.4,I5)
	WRITE(IO6,1104)(XMIN1(I),XMAX1(I),I=1,3),MDIR
 1104	FORMAT('0 CONTOUR MAP LIMIT AND DIRECTION',
	1 6F10.3,I8)
	READ(IO2, 104) (XMIN(I), XMAX(I), I = 1, 3), DX, NDIR
	IF(DX(1).GE.1.0.AND.DX(2).GE.1.0.AND.DX(3).GE.1.0)GO TO 1105
	GO TO 1106
 1105	DO 1107 I=1,3
	XMIN(I)=XMIN(I)/DX(I)
	XMAX(I)=XMAX(I)/DX(I)
 1107	DX(I)=1.0/DX(I)
 1106	CONTINUE
	WRITE(IO6, 105) (XMIN(I), XMAX(I), I = 1, 3),DX, NDIR
	IF(MDIR.EQ.NDIR)GO TO 1210
	GO TO (1201,1202,1203,1204,1205,1206),MDIR
 1207	WRITE(IO6,1209)MDIR,NDIR
 1209	FORMAT(1H ,I5,' & ',I5,' DIRECTIONS NOT COMPATABLE ')
	STOP
 1201	IF(NDIR.EQ.2)GO TO 1208
	GO TO 1207
 1202	IF(NDIR.EQ.1)GO TO 1208
	GO TO 1207
 1203	IF(NDIR.EQ.4)GO TO 1208
	GO TO 1207
 1204	IF(NDIR.EQ.3)GO TO 1208
	GO TO 1207
 1205	IF(NDIR.EQ.6)GO TO 1208
	GO TO 1207
 1206	IF(NDIR.NE.5)GO TO 1207
 1208	NDIR = MDIR
	MDIR = -MDIR
 1210	CONTINUE
 104	FORMAT(9F5.3, I8)
 105	FORMAT('0 INPUT MAP LIMITS AND DIRECTION '/' ',
	1 9F7.4, 5X, I3)
C      ROTATE PARAMETER LABELS TO STANDARD FORM
 3	CALL STRCON_ROT(XMIN(1),XMIN(2),XMIN(3),NDIR)
	CALL STRCON_ROT(XMAX(1),XMAX(2),XMAX(3),NDIR)
	CALL STRCON_ROT(DX(1), DX(2),DX(3),NDIR)
	CALL STRCON_ROT(XMIN1(1),XMIN1(2),XMIN1(3),NDIR)
	CALL STRCON_ROT(XMAX1(1),XMAX1(2),XMAX1(3),NDIR)
	CALL STRCON_ROT(A(1),A(2),A(3),NDIR)
	CALL STRCON_ROT(ANG(1), ANG(2),ANG(3),NDIR)
	CALL STRCON_ROT(XCHAR(1), XCHAR(2),XCHAR(3),NDIR)
	CALL STRCON_ROT(ACHAR(1), ACHAR(2),ACHAR(3),NDIR)
 4	CONTINUE
	IHEX = 0
	D1 = A(1)**2 + A(2)**2 + 2.0*A(1)*A(2)*COS(0.01745329*ANG(3))
	D2 = A(1)**2 + A(2)**2 - 2.0*A(1)*A(2)*COS(0.01745329*ANG(3))
	D1=D1/D2
	IF(D1.LT.1.0)D1=1.0/D1
	IF(D1.GT.1.7)IHEX=1
	DO 1220 I=1,3
	NN=IFIX(((XMIN1(I)-XMIN(I))/DX(I))+10.0)-10
	XMIN1(I)=XMIN(I)+NN*DX(I)
	NN=IFIX(((XMAX1(I)-XMIN(I))/DX(I))+10.99)-10
	XMAX1(I)=XMIN(I)+NN*DX(I)
 1220	CONTINUE
	WRITE(IO6,105)(XMIN1(I),XMAX1(I),I=1,3),(DX(I),I=1,3),NDIR
C      NOTE THAT ROTATING THE CHARACTER ARRAYS PERMITS
C      CORRECT LABELLING OF THE MAPS
C      CHECK LIMIT ORDER AND ARRAY SIZE
	NSDONE = 0
	DO 9 I = 1, 2
	IF(XMIN(I) .GE. XMAX(I)) GO TO 10
	IDIM(I) = (XMAX(I) - XMIN(I))/DX(I) + 1.01
 9	CONTINUE
	IF(XMIN(3) .LT. XMAX(3)) GO TO 11
 10	WRITE(IO6, 114)
 114	FORMAT('0 LIMIT AND DIRECTION CARD HAS A MAX VALUE LESS THAN',
	1 '  A MINIMUM, EXECUTION HALTING.')
	CALL STRCON_QUIT(NSDONE)
 11	CONTINUE
	WRITE(IO6, 201) IDIM
 201	FORMAT('0DIMENSIONS OF INPUT ARRAY ', I3, ' X ', I3)
	IF(IDIM(1)*IDIM(2) .LE. 10000) GO TO 25
	WRITE(IO6, 113)
 113	FORMAT('0MAP SIZE TOO LARGE--NOT ENOUGH SPACE IN ARRAYS MAP AND CO
	1ORD.')
	CALL STRCON_QUIT(NSDONE)
 25	CONTINUE
	READ(IO2, 106) RMIN, RMAX, RINC, NSEC, NEXT
	READ(IO2, 106) RMINN, RMAXN, RINCN
	MSEC=(XMAX1(3)-XMIN1(3))/DX(3) + 1.01
	IF(NSEC.GT.MSEC)NSEC=MSEC
 106	FORMAT(3F10.0, I5, I5)
	READ(IO2,126)SCALE,ISCALE,RSCX,RSCY,XFRAM,YFRAM,IFNFLD,
	1 XFNFLD,ISTRO,STRANG,MSEC
 126	FORMAT(F8.3,I5,4F8.3,I5,F8.4,I5,F8.4,I5)
	WRITE(IO6,371)SCALE,RSCX,RSCY,XFRAM,YFRAM,IFNFLD,
	1 XFNFLD,ISTRO,STRANG,MSEC
 371	FORMAT('  SCALE IS',F6.2,'MM/A.  ', /,
	1 '   DISTORTION CORRECTION IN VERT. AND HORIZ.=',2F10.3,/,
	2 '   VERT. AND HORIZ. FRAME =',2F10.3,' MM.',/,
	3 '   FAN FOLD PAPER?',I2,'  PAGE LENGTH =',F8.1,'MM.',/,
	4 '   STEREO ?',I2,' STEREO ANGLE =',F5.2,I5,' DENSITY SECTIONS/'
	5 'SECTION',//)
 388	CONTINUE
	IF( NDIR.EQ.1 .OR. NDIR.EQ.4 .OR. NDIR.EQ.5 ) STRANG = -STRANG
	MSEC2 = MSEC/2
	IF(MSEC2.EQ.0)MSEC2 = 1
	STRANG = DGTRD*STRANG*0.5
	STROX  = 0.00
	STROY  = 0.00
	SSTROX = 0.00
	SSTROY = 0.0
	IF( ISTRO.EQ.1 ) STROX = TAN(STRANG)
	IF( ISTRO.EQ.2 ) STROY = TAN(STRANG)
	IF(RINC  .EQ. 0) RINC = (RMAX - RMIN)/10
	NSEC1 = (XMAX(3)-XMIN(3))/DX(3)+1.01
	IF(NSEC.GT.NSEC1.OR.NSEC.EQ.0)NSEC = NSEC1
	NZMX = NSEC1
	NONEG = 0
	IF(RMINN.EQ.0.0.AND.RMAXN.EQ.0.0.AND.RINCN.EQ.0.0)NONEG = 1
	WRITE(IO6, 107) RMIN, RMAX, RINC, NSEC, NEXT, NMAXSD
	IF(NONEG.EQ.0)WRITE(IO6,108)RMINN,RMAXN,RINCN
	IF(NONEG.EQ.1)WRITE(IO6,109)
 109	FORMAT('0NO NEGATIVE CONTOURING CARRIED OUT')
 108	FORMAT('0WILL NEG CON FROM ',F6.0,' TO ',F6.0,' IN INCREMENTS OF '
	1 , F5.0)
 107	FORMAT('0WILL CONTOUR FROM ', F6.0, ' TO ', F6.0, ' IN INCREMENTS
	1 OF ', F5.0, ' FOR ', I2, ' SECTIONS.',/,1H ,' NEXT =',I5,/,1H ,
	2' NMAXSD = ',I5)
 40	CONTINUE
C
C     READ IN EXTENSION CARDS, AND CHECK THEM (IF NECESSARY)
C
	DO 501 I=1,3
	XMIN2(I)=XMIN(I)
	XMAX2(I)=XMAX(I)
 501	CONTINUE
	IF(NEXT .LE. 1) GO TO 520
	NDZ=0
	NEXT=0
	IF(XMAX(1)-XMIN(1)+DX(1).GT.0.999)GO TO 512
	WRITE(IO6,513)XCHAR(1),XMIN(1),XMAX(1),
	1 XMIN1(1),XMAX1(1),DX(1)
 513	FORMAT('0**** EXTENSION ALONG ',A1,' NOT POSSIBLE,',/,
	1 '  INPUT AND DESIRED LIMITS ARE',5F6.3)
	GO TO 516
 512	IF(XMIN1(1).LT.XMIN(1).OR.XMAX1(1).GT.XMAX(1))NEXT=1
 516	IF(XMAX(2)-XMIN(2)+DX(2).GT.0.999)GO TO 612
	WRITE(IO6,513)XCHAR(2),XMIN(2),XMAX(2),XMIN1(2),XMAX1(2),DX(2)
	GO TO 614
 612	IF(XMIN1(2).LT.XMIN(2).OR.XMAX1(2).GT.XMAX(2))NEXT=NEXT+2
 614	IF(XMIN1(3).GE.XMIN(3).AND.XMAX1(3).LE.XMAX(3))GO TO 510
	WRITE(IO6,511)XMIN(3),XMAX(3),XMIN1(3),XMAX1(3)
 511	FORMAT('0 E.D. SECTIONS',F8.2,' TO ',F8.2,' DO NOT ',
	1 'CONTAIN DESIRED SECTIONS',F8.2,' TO ',F8.2)
	STOP
 510	CONTINUE
	IINZ = IFIX( (XMIN1(3)-XMIN(3))/DX(3) + 0.1 )
	IINZ = IINZ - MSEC2
	IF( MSEC.EQ.1 ) IINZ = IINZ + 1
	IF( MDIR.LT.0 ) NDZ  = IDIM(2)
	IF( MDIR.GT.0 ) NDZ  = IDIM(1)
500 	CONTINUE
520 	CONTINUE
C      CHECK EXTENSION MODE, BE SURE FULL SECTION IS AVAILABLE
	MEXT = 1
	DO 5 I = 1, 2
 5	IF((1.0-(XMAX(I)-XMIN(I))).GE.(DX(I)+0.0001))MEXT=0
	IF(MEXT*NEXT .GT. 0) GO TO 7
	IF(MEXT .EQ. 0) GO TO 6
	WRITE(IO6, 110)
 110	FORMAT('0IT IS POSSIBLE TO EXTEND MAP TO FULL PLOTTING AREA BY USI
	1NG LATTICE TRANSLATIONS, HOWEVER THIS WILL NOT BE DONE.')
	GO TO 8
 6	WRITE(IO6, 111)
 111	FORMAT('0EXTENSION TO FULL PLOTTING AREA IMPOSSIBLE BECAUSE FULL S
	1ECTION IS NOT AVAILABLE.  EXECUTION CONTINUING WITHOUT EXTENSION')
	GO TO 8
 7	IF(NEXT .GE. 1) WRITE(IO6, 112)
 112	FORMAT('0MAP WILL BE EXTENDED TO FULL PLOTTING AREA BY LATTICE TRA
	1NSLATIONS.')
 8	CONTINUE
	IF(NEXT .GE. 1) NEXT = 3
	DO 530 I = 1, 3
	XMAX(I) = XMAX1(I)
	XMIN(I)=XMIN1(I)
530 	CONTINUE
C      SCALING SECTION
C      GET THE SHIFT FACTORS
C	CONVERT TO ANGSTROM SCALE
	DO 414 I=1,3
	DX(I)=DX(I)*A(I)
	XMIN(I)=XMIN(I)*A(I)
	XMAX(I)=XMAX(I)*A(I)
	XMIN1(I)=XMIN1(I)*A(I)
	XMAX1(I)=XMAX1(I)*A(I)
	XMIN2(I)=XMIN2(I)*A(I)
	XMAX2(I)=XMAX2(I)*A(I)
 414	CONTINUE
	SN1=SIN(ANG(1)*RC)
	CS1=COS(ANG(1)*RC)
	SN2=SIN(ANG(2)*RC)
	CS2=COS(ANG(2)*RC)
	SN3=SIN(ANG(3)*RC)
	CS3=COS(ANG(3)*RC)
	XSHX=1.0
	XSHY = CS3
	XSHZ = CS2
	YSHY=SN3
	YSHZ = (CS1-CS2*CS3)/SN2
	ZSHZ=SQRT((1.0-CS1**2-CS2**2-CS3**2+2.0*CS1*CS2*CS3)/SN3**2)
	WRITE(IO6,202) XSHX,XSHY,XSHZ,YSHY,YSHZ,ZSHZ
 202	FORMAT('0SHIFT FACTORS:  XSHX=', F6.3, ' XSHY=', F6.3,
	1' XSHZ=',F6.3,' YSHY=',F6.3,' YSHZ=',F6.3,' ZSHZ='F6.3)
C      GET MAXIMUM AND MINUMUM MAP COORDINATES
	XMPMIN = AMIN1(XSHY*XMIN(2), XSHY*XMAX(2), XSHZ*XMIN(3),
	1  XSHZ*XMAX(3), (XSHY*XMIN(2) + XSHZ*XMIN(3)),
	2  (XSHY*XMIN(2) + XSHZ*XMAX(3)), (XSHY*XMAX(2) + XSHZ*XMIN(3)),
	3  (XSHY*XMAX(2) + XSHZ*XMAX(3)))
	XMPMIN = XMPMIN + (XMIN(1) - DX(1))
	XMPMAX = AMAX1(XSHY*XMIN(2), XSHY*XMAX(2), XSHZ*XMIN(3),
	1XSHZ*XMAX(3), (XSHY*XMIN(2) + XSHZ*XMIN(3)),
	2  (XSHY*XMIN(2) + XSHZ*XMAX(3)), (XSHY*XMAX(2) + XSHZ*XMIN(3)),
	3  (XSHY*XMAX(2) + XSHZ*XMAX(3)))
	XMPMAX = XMPMAX + (XMAX(1) + DX(1))
	YMPMIN=AMIN1(YSHZ*XMIN(3),YSHZ*XMAX(3))+(XMIN(2)-DX(2))*YSHY
	YMPMAX=AMAX1(YSHZ*XMIN(3),YSHZ*XMAX(3))+(XMAX(2)+DX(2))*YSHY
	DO 300 IEX=1,3
	DELANG=ABS(ANG(IEX)-90.0)
	IF( DELANG.GT.0.01) GO TO 301
 300	CONTINUE
	XMPMAX=XMPMAX-DX(1)
	XMPMIN=XMPMIN+DX(1)
	YMPMAX=YMPMAX-DX(2)*YSHY
	YMPMIN=YMPMIN+DX(2)*YSHY
 301	CONTINUE
	WRITE(IO6, 115) XMPMIN, XMPMAX, YMPMIN, YMPMAX
 115	FORMAT('0A CARTESIAN COORDINATE SYSTEM FOR THE MAPS HAS BEEN ESTAB
	1LISHED'/' X RUNS FROM ', F7.2, ' TO ', F7.2/' Y RUNS FROM ',
	2  F7.2, ' TO ', F7.2/' THE COORDINATES OF EACH GRID POINT W'
	3 'ILL BE TRANSFORMED TO THIS SYSTEM BEFORE CONTOURING.')
C      CHECK TO SEE THAT CONTOURING ARRAY DIMESNIONS AREN'T EXCEEDED
	JDIM(1) = (XMPMAX - XMPMIN)/DX(1) + 1.01
	JDIM(2) = (YMPMAX - YMPMIN)/(DX(2)*YSHY) + 1.01
	WRITE(IO6, 203) JDIM
 203	FORMAT('0DIMENSIONS OF CONTOURING ARRAY ', I3, ' X ', I3)
	IF((JDIM(1).LE.109).AND.(JDIM(2).LE.109)) GO TO 26
	WRITE(IO6, 116)
 116	FORMAT('0 NOT ENOUGH SPACE IN CONTOURING ARRAY, 109 X 109 MAX.')
	CALL STRCON_QUIT(NSDONE)
 26	CONTINUE
C      GET SCALE FACTOR
	WRITE(IO6, 117) ISCALE
 117	FORMAT('0GRID SPACING IS ', I3, 'A.')
	CALL STRCON_MODEL1(NDIR)
	CLOSE(UNIT=IO2)
	RETURN
	END

	SUBROUTINE STRCON_VDWDRW(SSTROX,SSTROY,STML0,STMT0,MSEC,
     +  STROSN,SCX,SCY)

	COMMON /IOSET/ IO1,IO2,IO3,IO4,IO5,IO6,IO7,IO8
	COMMON /MOD/ AM(3,3,36),AT(3,36),XMIN3(3),XMAX3(3),
	1 DZ3,DZ3N,DZ2,DZ2N,DMX(4),CHAR(4,10),NSYM,PWIDE(4),ICHR(4),
	2 HT(4),HT2(4),SCALE,RAD(4,10),RDEL(4),IMASK(36,4),
	3 NSET(4,20),ATMFIL(36,4),VRAD(10),VTYP(10)
	COMMON /DAT/ TI(10), XMIN(3), XMAX(3), DX(3), XCHAR(3), ACHAR(3),
	1  IDIM(2),JDIM(2),NEXT,XSHX,XSHY,XSHZ,YSHY,YSHZ,ZSHZ,XMPMIN,
	2  XMPMAX, YMPMIN, YMPMAX, ISCALE, NSDONE, NIP, ZMAP, A,
	3  RMIN, RMAX, RINC, NSEC, ANG, RC
	DIMENSION XI(3,2500),KI(2500),X(3),X1(3)
	DIMENSION YIN(4,2),YTR(4,2)
	LOGICAL*1 KI,ATMFIL
	INTEGER*4 PWIDE,CHAR,VTYP
	IF(NSYM.EQ.0)RETURN
	CALL STRCON_CIRCLE
	DO 100 II=1,4
	IF(DMX(II).EQ.0.0.AND.HT(II).EQ.0.0)GO TO 100
	DMX2=DMX(II)**2
	XMIN3(3)=-DMX(II)-DZ2
	XMAX3(3)=DMX(II)+DZ2
	ZMAPO=ZMAP*ZSHZ
	J=0
	OPEN(UNIT=IO4,FILE=ATMFIL(1,II),TYPE='OLD',FORM='FORMATTED')
 10	CONTINUE
	READ(IO4,11,END=99)N,IK,X(1),X(2),X(3)
 11	FORMAT(I5,10X,I3,3F10.4)
	DO 110 K=1,19,2
	IF((N.GE.NSET(II,K)).AND.(N.LE.NSET(II,K+1)))GO TO 120
 110	CONTINUE
	GO TO 10
 120	CONTINUE
	IF(IK.LE.0)GO TO 10
	DO 20 I=1,NSYM
	IF(IMASK(I,II).EQ.0)GO TO 20
	DO 22 I1=1,3
	X1(I1)=AT(I1,I)
	DO 24 I2=1,3
	X1(I1)=X1(I1)+AM(I1,I2,I)*X(I2)
 24	CONTINUE
 22	CONTINUE
	X1(3)=X1(3)-ZMAPO
	DO 30 I1=1,3
	IF(X1(I1).LT.XMIN3(I1).OR.X1(I1).GT.XMAX3(I1))GO TO 20
 30	CONTINUE
	J=J+1
	IF(J.GT.2500)GO TO 38
	DO 32 I1=1,3
	XI(I1,J)=X1(I1)
 32	CONTINUE
	IF(IK.GT.10)IK=10
	KI(J)=IK
 20	CONTINUE
	GO TO 10
 38	WRITE(IO6,39)
 39	FORMAT('0*** TOO MANY ATOMS ***')
	CALL STRCON_QUIT(0)
 99	CONTINUE
	CLOSE(UNIT=IO4)
	JN=J
C
C	JN ATOMS TO BE INCLUDED HAVE BEEN COLLECTED.
C
	IF(JN.EQ.0)GO TO 100
	DO 1000 ISEC=1,MSEC
	DZ=0.5*FLOAT(2*ISEC-1-MSEC)
	STML=STML0-DZ*SSTROY*STROSN
	STMT=STMT0-DZ*SSTROX*STROSN
	DZ=DZ*DX(3)*ZSHZ
	CALL STRCON_OFFSET(STMT,SCX,STML,SCY)
	DO 40 J=1,JN
C
C	ESTABLISH LIMITS FOR ATOM J.
C	NOTE:	CIRCLE IS DEVIDED INTO TO HALVES TOP AND BOTTOM
C		FOR TOP   - YIN IS LEFT  AND YTR IS RIGHT
C		FOR BOTTOM- YIN IS RIGHT AND YTR IS LEFT
C
	X0=XI(1,J)
	Y0=XI(2,J)
	IK=KI(J)
	IF(VTYP(IK).EQ.0)GO TO 40
	R0=VRAD(IK)**2-(XI(3,J)-DZ)**2
	IF(R0.LE.0.0)GO TO 40
	R0=SQRT(R0)
	XIN=X0-R0
	YIN(1,1)=Y0-R0
	XTR=X0+R0
	YTR(1,1)=Y0+R0
	YIN(1,2)=Y0+R0
	YTR(1,2)=Y0-R0
	YIN(2,1)=0.0
	YTR(2,1)=0.0
	YIN(3,1)=0.0
	YTR(3,1)=0.0
	YIN(2,2)=0.0
	YTR(2,2)=0.0
	YIN(3,2)=0.0
	YTR(3,2)=0.0
	YIN(4,1)=0.0
	YTR(4,1)=0.0
	YIN(4,2)=0.0
	YTR(4,2)=0.0
 60	IF(YIN(1,1).GT.YMPMIN)GO TO 62
	IF(YTR(1,1).LE.YMPMIN)GO TO 40
	YIN(1,1)=YMPMIN
	YTR(1,2)=YMPMIN
	IF(YMPMIN.LT.Y0)GO TO 62
	DDX=SQRT(R0**2-(YIN(1,1)-Y0)**2)
	XIN=X0-DDX
	XTR=X0+DDX
 62	IF(YTR(1,1).LT.YMPMAX)GO TO 64
	IF(YIN(1,1).GT.YMPMAX)GO TO 40
	YTR(1,1)=YMPMAX
	YIN(1,2)=YMPMAX
	IF(Y0.LE.YMPMAX)GO TO 64
	DDX=SQRT(R0**2-(YTR(1,1)-Y0)**2)
	XIN=X0-DDX
	XTR=X0+DDX
 64	IF(XIN.GT.XMPMIN)GO TO 66
	IF(XTR.LT.XMPMIN)GO TO 40
	XIN=XMPMIN
	IF(X0.LT.XMPMIN)GO TO 65
	YTR(2,1)=YTR(1,1)
	DY=SQRT(R0**2-(X0-XMPMIN)**2)
	YTR(1,1)=Y0-DY
	YIN(2,1)=Y0+DY
	GO TO 66
 65	YIN(1,1)=0.0
	YTR(1,1)=0.0
	DY=SQRT(R0**2-(X0-XMPMIN)**2)
	YTR(1,2)=Y0-DY
	YIN(1,2)=Y0+DY
 66	IF(XTR.LT.XMPMAX)GO TO 68
	IF(XIN.GT.XMPMAX)GO TO 40
	XTR=XMPMAX
	IF(X0.GT.XMPMAX)GO TO 67
	YTR(2,2)=YTR(1,2)
	DY=SQRT(R0**2-(XMPMAX-X0)**2)
	YIN(2,2)=Y0-DY
	YTR(1,2)=Y0+DY
	GO TO 68
 67	YTR(1,2)=0.0
	YIN(1,2)=0.0
	DY=SQRT(R0**2-(XMPMAX-X0)**2)
	YIN(1,1)=Y0-DY
	YTR(1,1)=Y0+DY
 68	CONTINUE
C
C	LIMITS ESTABLISHED.
C	DETERMINE BONDED ATOMS.
C
	DO  160 J2=1,JN
	IF(J.EQ.J2)GO TO 160
	D=0.0
	DO 162 K=1,3
	D=D+(XI(K,J)-XI(K,J2))**2
	IF(D.GT.DMX2)GO TO 160
 162	CONTINUE
	D=SQRT(D)
	KJ1=KI(J)
	KJ2=KI(J2)
	IF(D.GT.RDEL(II)+ABS(RAD(II,KJ1))+ABS(RAD(II,KJ2)))GO TO 160
C
C	ATOMS BONDED.
C
	R2=VRAD(KJ2)**2-(XI(3,J2)-DZ)**2
	IF(R2.LE.0.0)GO TO 160
	R2=SQRT(R2)
	X2=XI(1,J2)
	Y2=XI(2,J2)
	RR=(X2-X0)**2+(Y2-Y0)**2
	RR=SQRT(RR)
C
C	IF ATOM-0 TOTALLY INSIDE ATOM-2 GO TO NEXT ATOM.
C
	IF((R0+RR).LE.R2)GO TO 40
C
C	FIND INTERSECTION OF CIRCLES. NOTE X3 < OR = X4.
C
	CALL STRCON_INTRSC(X0,Y0,R0,X2,Y2,R2,X3,Y3,X4,Y4)
	IF(X3.EQ.0.0.AND.X4.EQ.0.0.AND.Y3.EQ.0.0.AND.Y4.EQ.0.0)GO TO 160
C
C	IF ATOM-0 INSIDE CIRCLE OF ATOM-2 GO TO 240.
C
	IF(RR.LT.R2)GO TO 240
C
C	ATOM-0 OUTSIDE CIRCLE OF ATOM-2.
C	IF ATOM-2 RIGHT OF ATOM-0 GO TO 170.
	IF(Y2.GT.Y0)GO TO 170
	IF(X4.LT.X0)GO TO 180
	IF(X3.GT.X0)GO TO 190
C
C	ATOM-2 LEFT OF ATOM-0.
C
	DO 202 K=1,4
	IF(YIN(K,1).EQ.0.0.AND.YTR(K,1).EQ.0.0)GO TO 204
	IF(YIN(K,1).GE.Y3)GO TO 204
	YIN(K,1)=Y3
	IF(YIN(K,1).LT.YTR(K,1))GO TO 204
C
C	YIN AND YTR HAVE CROSSED, THIS PORTION BLANKED OUT.
C
	YIN(K,1)=0.0
	YTR(K,1)=0.0
 204	IF(YIN(K,2).EQ.0.0.AND.YTR(K,2).EQ.0.0)GO TO 202
	IF(YTR(K,2).GE.Y4)GO TO 202
	YTR(K,2)=Y4
	IF(YTR(K,2).LT.YIN(K,2))GO TO 202
C
C	YIN AND YTR HAVE CROSSED OVER, THIS PORTION OBLITERATED.
C
	YIN(K,2)=0.0
	YTR(K,2)=0.0
 202	CONTINUE
	GO TO 164
 180	CONTINUE
C
C	ATOM-2 ABOVE ATOM-0.
C	PLACE Y3 ON LEFT AND Y4 ON RIGHT.
C
	Y=Y3
	Y3=Y4
	Y4=Y
	Y=X3
	X3=X4
	X4=Y
	GO TO 210
 190	CONTINUE
C
C	ATOM-2 BELOW ATOM-0 (>X).
C	PLACE Y4 ON LEFT AND Y3 ON RIGHT.
C
	Y=Y3
	Y3=Y4
	Y4=Y
	Y=X3
	X3=X4
	X4=Y
	GO TO 220
 170	CONTINUE
C
C	ATOM-2 RIGHT OF ATOM-0.
C
	IF(X4.LT.X0)GO TO 210
	IF(X3.GT.X0)GO TO 220
	DO 232 K=1,4
	IF(YIN(K,1).EQ.0.0.AND.YTR(K,1).EQ.0.0)GO TO 234
	IF(YTR(K,1).LE.Y3)GO TO 234
	YTR(K,1)=Y3
	IF(YTR(K,1).GT.YIN(K,1))GO TO 234
C
C	YIN AND YTR CROSSED OVER, PORTION OBLITERATED.
C
	YIN(K,1)=0.0
	YTR(K,1)=0.0
 234	IF(YIN(K,2).EQ.0.0.AND.YTR(K,2).EQ.0.0)GO TO 232
	IF(YIN(K,2).LE.Y4)GO TO 232
	YIN(K,2)=Y4
	IF(YTR(K,2).LT.YIN(K,2))GO TO 232
C
C	YIN AND YTR CROSSED OVER, PORTION OBLITERATED.
C
	YIN(K,2)=0.0
	YTR(K,2)=0.0
 232	CONTINUE
	GO TO 164
 210	CONTINUE
C
C	ATOM-2 ABOVE ATOM-0 (X2&X3&X4 < X0)
C	ALL IN TOP HALF - Y3 LEFT, Y4 RIGHT
	DO 212 K=1,4
	IF(YIN(K,1).EQ.0.0.AND.YTR(K,1).EQ.0.0)GO TO 212
	IF(YIN(K,1).LT.Y3.AND.YTR(K,1).GT.Y4)GO TO 214
	IF(Y3.LT.YTR(K,1).AND.Y4.GT.YTR(K,1))YTR(K,1)=Y3
	IF(Y3.LT.YIN(K,1).AND.Y4.GT.YIN(K,1))YIN(K,1)=Y4
	IF(YIN(K,1).LT.YTR(K,1))GO TO 212
C
C	YIN AND YTR CROSSED OVER, PORTION OBLITERATED.
C
	YIN(K,1)=0.0
	YTR(K,1)=0.0
	GO TO 212
C
C	SEGMENT MUST BE DEVIDED.
C
 214	K1=K
 216	K1=K1+1
	IF(K1.GT.4)GO TO 212
	IF(YIN(K1,1).NE.0.0.OR.YTR(K1,1).NE.0.0)GO TO 216
	YTR(K1,1)=YTR(K,1)
	YTR(K,1)=Y3
	YIN(K1,1)=Y4
 212	CONTINUE
	GO TO 164
 220	CONTINUE
C
C	ALL IN BOTTOM HALF - Y4 LEFT, Y3 RIGHT
C
	DO 222 K=1,4
	IF(YIN(K,2).EQ.0.0.AND.YTR(K,2).EQ.0.0)GO TO 222
	IF(YIN(K,2).GT.Y3.AND.YTR(K,2).LT.Y4)GO TO 224
	IF(Y3.GT.YIN(K,2).AND.Y4.LT.YIN(K,2))YIN(K,2)=Y4
	IF(Y3.GT.YTR(K,2).AND.Y4.LT.YTR(K,2))YTR(K,2)=Y3
	IF(YTR(K,2).LT.YIN(K,2))GO TO 222
C
C	YIN AND YTR CROSSED OVER, PORTION OBLITERATED.
C
	YTR(K,2)=0.0
	YIN(K,2)=0.0
	GO TO 222
C
C	SEGMENT MUST BE DEVIDED.
C
 224	K1=K
 226	K1=K1+1
	IF(K1.GT.4)GO TO 222
	IF(YIN(K1,2).NE.0.0.OR.YTR(K1,2).NE.0.0)GO TO 226
	YTR(K1,2)=YTR(K,2)
	YTR(K,2)=Y3
	YIN(K1,2)=Y4
 222	CONTINUE
	GO TO 164
C
C	ATOM-0 INSIDE CIRCLE OF ATOM-2.
C
 240	IF(Y2.LT.Y0)GO TO 270
C
C	ATOM-2 RIGHT OF ATOM-0.
C
	IF(X4.LT.X0)GO TO 280
	IF(X3.GT.X0)GO TO 290
C
C	ATOM-2 RIGHT OF ATOM-0.
C
	DO 302 K=1,4
	IF(YIN(K,1).EQ.0.0.AND.YTR(K,1).EQ.0.0)GO TO 304
	IF(YTR(K,1).LE.Y3)GO TO 304
	YTR(K,1)=Y3
	IF(YIN(K,1).LT.YTR(K,1))GO TO 304
C
C	YIN AND YTR CROSSED OVER, PORTION OBLITERATED.
C
	YIN(K,1)=0.0
	YTR(K,1)=0.0
 304	IF(YIN(K,2).EQ.0.0.AND.YTR(K,2).EQ.0.0)GO TO 302
	IF(YIN(K,2).LE.Y4)GO TO 302
	YIN(K,2)=Y4
	IF(YTR(K,2).LT.YIN(K,2))GO TO 302
C
C	YIN AND YTR CROSSED OVER, PORTION OBLITERATED.
C
	YIN(K,2)=0.0
	YTR(K,2)=0.0
 302	CONTINUE
	GO TO 164
 280	CONTINUE
C
C	ATOM-2 ABOVE ATOM-0 (X2&X3&X4 < X0)
C	PLACE Y3 ON LEFT AND Y4 ON RIGHT.
C
	Y=Y3
	Y3=Y4
	Y4=Y
	Y=X3
	X3=X4
	X4=Y
	GO TO 310
 290	CONTINUE
C
C	INTERSECTION BELOW ATOM-0 ( X3&X4 > X0 )
C	PLACE Y4 ON LEFT AND Y3 ON RIGHT.
C
	Y=Y3
	Y3=Y4
	Y4=Y
	Y=X3
	X3=X4
	X4=Y
	GO TO 320
 270	CONTINUE
C
C	ATOM-2 LEFT OF ATOM-0.
C
	IF(X4.LT.X0)GO TO 310
	IF(X3.GT.X0)GO TO 320
C
C	ATOM-2 LEFT OF ATOM-0.
C
	DO 332 K=1,4
	IF(YIN(K,1).EQ.0.0.AND.YTR(K,1).EQ.0.0)GO TO 334
	IF(YIN(K,1).GE.Y3)GO TO 334
	YIN(K,1)=Y3
	IF(YTR(K,1).GT.YIN(K,1))GO TO 334
C
C	YIN AND YTR CROSSED OVER, PORTION OBLITERATED.
C
	YIN(K,1)=0.0
	YTR(K,1)=0.0
 334	IF(YIN(K,2).EQ.0.0.AND.YTR(K,2).EQ.0.0)GO TO 332
	IF(YTR(K,2).GE.Y4)GO TO 332
	YTR(K,2)=Y4
	IF(YTR(K,2).LT.YIN(K,2))GO TO 332
C
C	YIN AND YTR CROSSED OVER, PORTION OBLITERATED.
C
	YIN(K,2)=0.0
	YTR(K,2)=0.0
 332	CONTINUE
	GO TO 164
 310	CONTINUE
C
C	ALL IN TOP HALF - Y3 LEFT, Y4 RIGHT
C
	IF(X2.LT.X0)GO TO 180
C
C	BOTTOM HALF OBLITERATED.
C
	DO 312 K=1,4
	YIN(K,2)=0.0
	YTR(K,2)=0.0
	IF(YIN(K,1).EQ.0.0.AND.YTR(K,1).EQ.0.0)GO TO 312
	IF(YIN(K,1).GT.Y3.AND.YTR(K,1).LT.Y4)GO TO 312
	IF(Y3.GT.YIN(K,1).AND.Y3.LT.YTR(K,1))YIN(K,1)=Y3
	IF(Y4.GT.YIN(K,1).AND.Y4.LT.YTR(K,1))YTR(K,1)=Y4
	IF(YIN(K,1).LT.YTR(K,1))GO TO 312
C
C	YIN AND YTR CROSSED OVER, PORTION OBLITERATED.
C
	YIN(K,1)=0.0
	YTR(K,1)=0.0
 312	CONTINUE
	GO TO 164
 320	CONTINUE
C
C	ALL IN BOTTOM HALF - Y3 LEFT, Y4 RIGHT
C
	IF(X2.GT.X0)GO TO 190
	DO 322 K=1,4
C
C	TOP HALF OBLITERATED.
C
	YIN(K,1)=0.0
	YTR(K,1)=0.0
	IF(YIN(K,2).EQ.0.0.AND.YTR(K,2).EQ.0.0)GO TO 322
	IF(YIN(K,2).LT.Y4.AND.YTR(K,2).GT.Y3)GO TO 322
	IF(Y3.LT.YIN(K,2).AND.Y3.GT.YTR(K,2))YIN(K,2)=Y3
	IF(Y4.LT.YIN(K,2).AND.Y4.GT.YTR(K,2))YTR(K,2)=Y4
	IF(YTR(K,2).LT.YIN(K,2))GO TO 322
C
C	YIN AND YTR CROSSED OVER, PORTION OBLITERATED.
C
	YTR(K,2)=0.0
	YIN(K,2)=0.0
 322	CONTINUE
 164	CONTINUE
	DO 350 K=1,3
	IF(YTR(K,1).GT.YIN(K,1))GO TO 354
	YIN(K,1)=YIN(K+1,1)
	YIN(K+1,1)=0.0
	YTR(K,1)=YTR(K+1,1)
	YTR(K+1,1)=0.0
 354	IF(YIN(K,2).GT.YTR(K,2))GO TO 350
	YIN(K,2)=YIN(K+1,2)
	YIN(K+1,2)=0.0
	YTR(K,2)=YTR(K+1,2)
	YTR(K+1,2)=0.0
 350	CONTINUE
 160	CONTINUE
	DO 90 K=1,4
	Y1=YIN(K,1)
	Y2=YTR(K,1)
	Y3=YIN(K,2)
	Y4=YTR(K,2)
	IF(Y1.LT.YMPMIN)Y1=YMPMIN
	IF(Y2.GT.YMPMAX)Y2=YMPMAX
	IF(Y3.GT.YMPMAX)Y3=YMPMAX
	IF(Y4.LT.YMPMIN)Y4=YMPMIN
	IF(YIN(K,1).NE.0.0)CALL STRCON_CIRDRW(X0,Y0,R0,Y1,Y2,VTYP(IK))
	IF(YIN(K,2).NE.0.0)CALL STRCON_CIRDRW(X0,Y0,R0,Y3,Y4,VTYP(IK))
	YIN(K,1)=0.0
	YIN(K,2)=0.0
	YTR(K,1)=0.0
	YTR(K,2)=0.0
 90	CONTINUE
 40	CONTINUE
 1000	CONTINUE
 100	CONTINUE
	RETURN
	END

	SUBROUTINE STRCON_CIRCLE
	DIMENSION XUN(21),YUN(21)
	PI=3.141592654
	DT=PI/20.
	DO 20 I=1,21
	YUN(I)=-COS(DT*FLOAT(I-1))
	XUN(I)=SQRT(1.0-YUN(I)**2)
 20	CONTINUE
	RETURN
C
C
	ENTRY STRCON_CIRDRW(X0,Y0,R0,Y1,Y2,ITYP)
	IF(Y2.LT.Y1)GO TO 10
	T1=((Y1-Y0)/R0)
	T2=((Y2-Y0)/R0)
	IF(T1.LT.-1.0)T1=-1.0
	IF(T1.GT.1.0)T1=1.0
	IF(T2.LT.-1.0)T2=-1.0
	IF(T2.GT.1.0)T2=1.0
	T1=ACOS(-T1)
	T2=ACOS(-T2)
	I1=(T1/DT)+2.0
	I2=(T2/DT)+1.0
	I3=(0.5+1.5/R0)
	IF(I3.LT.1)I3=1
	IF(I3.GT.4)I3=4
	Y=Y1
	X=(R0**2-(Y-Y0)**2)
	IF(X.LT.0.0)X=0.0
	X=X0-SQRT(X)
	CALL STRCON_QPLOT(X,Y,+13)
	GO TO (12,112,212,312)ITYP
 12	CONTINUE
	IF(I2.LE.I1)GO TO 40
	DO 30 I=I1,I2,I3
	X=X0-R0*XUN(I)
	Y=Y0+R0*YUN(I)
	CALL STRCON_QPLOT(X,Y,12)
 30	CONTINUE
 40	CONTINUE
	Y=Y2
	X=(R0**2-(Y-Y0)**2)
	IF(X.LT.0.0)X=0.0
	X=X0-SQRT(X)
	CALL STRCON_QPLOT(X,Y,12)
	RETURN
C
C
 112	CONTINUE
	IF(I2.LE.I1)GO TO 140
	DO 130 I=I1,I2,I3
	CALL STRCON_QPLOT(X,Y,12)
	X=X0-R0*XUN(I)
	Y=Y0+R0*YUN(I)
	CALL STRCON_QPLOT(X,Y,13)
 130	CONTINUE
 140	CONTINUE
	Y=Y2
	X=(R0**2-(Y-Y0)**2)
	IF(X.LT.0.0)X=0.0
	X=X0-SQRT(X)
	CALL STRCON_QPLOT(X,Y,13)
	CALL STRCON_QPLOT(X,Y,12)
	RETURN
C
C
 212	CONTINUE
	IF(I2.LE.I1)GO TO 240
	I4=13
	I5=12
	DO 230 I=I1,I2,I3
	X=X0-R0*XUN(I)
	Y=Y0+R0*YUN(I)
	CALL STRCON_QPLOT(X,Y,I4)
	I6=I4
	I4=I5
	I5=I6
 230	CONTINUE
 240	CONTINUE
	Y=Y2
	X=(R0**2-(Y-Y0)**2)
	IF(X.LT.0.0)X=0.0
	X=X0-SQRT(X)
	CALL STRCON_QPLOT(X,Y,12)
	RETURN
C
C
 312	CONTINUE
	K=-1
	IF(I2.LE.I1)GO TO 340
	DO 330 I=I1,I2,I3
	K=-1*K
	X=X0-(R0+K*0.1)*XUN(I)
	Y=Y0+(R0+K*0.1)*YUN(I)
	CALL STRCON_QPLOT(X,Y,12)
 330	CONTINUE
 340	CONTINUE
	Y=Y2
	X=(R0**2-(Y-Y0)**2)
	IF(X.LT.0.0)X=0.0
	X=X0-SQRT(X)
	CALL STRCON_QPLOT(X,Y,12)
	RETURN
C
C
 10	CONTINUE
	T2=((Y1-Y0)/R0)
	T1=((Y2-Y0)/R0)
	IF(T1.LT.-1.0)T1=-1.0
	IF(T1.GT.1.0)T1=1.0
	IF(T2.LT.-1.0)T2=-1.0
	IF(T2.GT.1.0)T2=1.0
	T1=ACOS(-T1)
	T2=ACOS(-T2)
	I1=(T1/DT)+2.0
	I2=(T2/DT)+1.0
	I3=(0.5+1.5/R0)
	IF(I3.LT.1)I3=1
	IF(I3.GT.4)I3=4
	Y=Y2
	X=(R0**2-(Y-Y0)**2)
	IF(X.LT.0.0)X=0.0
	X=X0+SQRT(X)
	CALL STRCON_QPLOT(X,Y,+13)
	GO TO (14,114,214,314)ITYP
 14	CONTINUE
	IF(I2.LE.I1)GO TO 60
	DO 50 I=I1,I2,I3
	X=X0+R0*XUN(I)
	Y=Y0+R0*YUN(I)
	CALL STRCON_QPLOT(X,Y,12)
 50	CONTINUE
 60	CONTINUE
	Y=Y1
	X=(R0**2-(Y-Y0)**2)
	IF(X.LT.0.0)X=0.0
	X=X0+SQRT(X)
	CALL STRCON_QPLOT(X,Y,12)
	RETURN
C
C
 114	CONTINUE
	IF(I2.LE.I1)GO TO 160
	DO 150 I=I1,I2,I3
	CALL STRCON_QPLOT(X,Y,12)
	X=X0+R0*XUN(I)
	Y=Y0+R0*YUN(I)
	CALL STRCON_QPLOT(X,Y,13)
 150	CONTINUE
 160	CONTINUE
	Y=Y1
	X=(R0**2-(Y-Y0)**2)
	IF(X.LT.0.0)X=0.0
	X=X0+SQRT(X)
	CALL STRCON_QPLOT(X,Y,13)
	CALL STRCON_QPLOT(X,Y,12)
	RETURN
C
C
 214	CONTINUE
	IF(I2.LE.I1)GO TO 260
	I4=13
	I5=12
	DO 250 I=I1,I2,I3
	X=X0+R0*XUN(I)
	Y=Y0+R0*YUN(I)
	CALL STRCON_QPLOT(X,Y,I4)
	I6=I4
	I4=I5
	I5=I6
 250	CONTINUE
 260	CONTINUE
	Y=Y1
	X=(R0**2-(Y-Y0)**2)
	IF(X.LT.0.0)X=0.0
	X=X0+SQRT(X)
	CALL STRCON_QPLOT(X,Y,12)
	RETURN
C
C
 314	CONTINUE
	K=-1
	IF(I2.LE.I1)GO TO 360
	DO 350 I=I1,I2,I3
	K=-1*K
	X=X0+(R0+K*0.1)*XUN(I)
	Y=Y0+(R0+K*0.1)*YUN(I)
	CALL STRCON_QPLOT(X,Y,12)
 350	CONTINUE
 360	CONTINUE
	Y=Y1
	X=(R0**2-(Y-Y0)**2)
	IF(X.LT.0.0)X=0.0
	X=X0+SQRT(X)
	CALL STRCON_QPLOT(X,Y,12)
	RETURN
	END

	SUBROUTINE STRCON_INTRSC(X1,Y1,R1,X2,Y2,R2,X3,Y3,X4,Y4)
	IF(Y1.EQ.Y2)GO TO 30
	A2=0.5*(R1**2-X1**2-Y1**2-R2**2+X2**2+Y2**2)
	A2=A2/(Y1-Y2)
	A1=-(X1-X2)/(Y1-Y2)
	A=1+A1**2
	B=-2.0*(X1+A1*(A2+Y1))
	C=X1**2+(A2+Y1)**2-R1**2
	D=B**2-4.0*A*C
	IF(D.LT.0.0)GO TO 10
	D=SQRT(D)
	X3=(-B-D)/(2.0*A)
	X4=(-B+D)/(2.0*A)
	Y3=X3*A1-A2
	Y4=X4*A1-A2
	RETURN
 10	X3=0.0
	Y3=0.0
	X4=0.0
	Y4=0.0
	RETURN
 30	CONTINUE
	X3=R2**2-R1**2-X2**2+X1**2
	X3=0.5*X3/(X1-X2)
	X4=X3
	A=(R1**2-(X3-X1)**2)
	IF(A.LT.0.0)GO TO 10
	A=SQRT(A)
	Y3=Y1-A
	Y4=Y1+A
	RETURN
	END

	SUBROUTINE STRCON_MODEL(STROX,STROY)
C
C	Subroutine to mark atom positions and draw bonds
C
C
	COMMON /IOSET/ IO1,IO2,IO3,IO4,IO5,IO6,IO7,IO8
C
	COMMON /MOD/ AM(3,3,36),AT(3,36),XMIN3(3),XMAX3(3),
	1 DZ3,DZ3N,DZ2,DZ2N,DMX(4),CHAR(4,10),NSYM,PWIDE(4),ICHR(4),
	2 HT(4),HT2(4),SCALE,RAD(4,10),RDEL(4),IMASK(36,4),
	3 NSET(4,20),ATMFIL(36,4)
C
	COMMON /DAT/ TI(10), XMIN(3), XMAX(3), DX(3), XCHAR(3), ACHAR(3),
	1  IDIM(2),JDIM(2),NEXT,XSHX,XSHY,XSHZ,YSHY,YSHZ,ZSHZ,XMPMIN,
	2  XMPMAX, YMPMIN, YMPMAX, ISCALE, NSDONE, NIP, ZMAP, A,
	3  RMIN, RMAX, RINC, NSEC, ANG, RC
C
	DIMENSION XI(3,2500),KI(2500),X(3),X1(3)
	LOGICAL*1 KI,ATMFIL
	INTEGER*4 PWIDE,CHAR, LABLL
     	CHARACTER LABEL*4
     	EQUIVALENCE( LABEL, LABLL )
C
	IF(NSYM.EQ.0)RETURN
C
C	master loop over seperate atom sets.
C
	DO 100 II=1,4
	IF(DMX(II).EQ.0.0.AND.HT(II).EQ.0.0)GO TO 100
	DMX2=DMX(II)**2
	XMIN3(3)=-DMX(II)-DZ2
	XMAX3(3)=DMX(II)+DZ2
	ZMAPO=ZMAP*ZSHZ
	J=0
	OPEN(UNIT=IO4,FILE=ATMFIL(1,II),TYPE='OLD',FORM='FORMATTED')
C
 10	CONTINUE
	READ(IO4,11,END=99)N,IK,X(1),X(2),X(3)
 11	FORMAT(I5,10X,I3,3F10.4)
	DO K=1,19,2
	   IF((N.GE.NSET(II,K)).AND.(N.LE.NSET(II,K+1)))GO TO 120
	END DO
	GO TO 10
 120	CONTINUE
	IF(IK.LE.0)GO TO 10
	DO 20 I=1,NSYM
	   IF(IMASK(I,II).EQ.0)GO TO 20
	   DO I1=1,3
	      X1(I1)=AT(I1,I)
	      DO I2=1,3
	         X1(I1)=X1(I1)+AM(I1,I2,I)*X(I2)
	      END DO
	   END DO
	   X1(3)=X1(3)-ZMAPO
	   DO I1=1,3
	      IF(X1(I1).LT.XMIN3(I1).OR.X1(I1).GT.XMAX3(I1))GO TO 20
	   END DO
	   J=J+1
	   IF(J.GT.2500)GO TO 38
	   DO I1=1,3
	      XI(I1,J)=X1(I1)
	   END DO
	   IF(IK.GT.10)IK=10
	   KI(J)=IK
c	   write( 85, 11) j, ki(j), xi(1,j), xi(2,j), xi(3,j)
 20	CONTINUE
	GO TO 10
C
 38	WRITE(IO5,39)
 39	FORMAT('0*** TOO MANY ATOMS ***')
	CALL STRCON_QUIT(0)
 99	CONTINUE
	CLOSE(UNIT=IO4)
C
C	mark atom position.
C
	JN=J
	IF(JN.EQ.0)GO TO 100
C
	XMPMIN1 = XMPMIN + 1.2 * HT2(II) / SCALE
	XMPMAX1 = XMPMAX - 1.2 * HT2(II) / SCALE
	YMPMIN1 = YMPMIN + HT2(II) / SCALE
	YMPMAX1 = YMPMAX - HT2(II) / SCALE
	DO J=1,JN
	   IF( XI(3,J).GE.DZ2N .AND. XI(3,J).LE.DZ2 .AND.
	1  XI(1,J).GE.XMPMIN1 .AND. XI(1,J).LE.XMPMAX1 .AND.
	2  XI(2,J).GE.YMPMIN1 .AND. XI(2,J).LE.YMPMAX1 .AND.
	3  ICHR(II).NE.0 ) THEN
	      IK=KI(J)
     	      IF( CHAR(II,IK).NE.0 ) THEN
	      XC=XI(1,J)+STROX*XI(3,J)
	      YC=XI(2,J)+STROY*XI(3,J)
     	      LABLL = CHAR(II,IK)
     	      CALL STRCON_TEXT(XC,YC,LABEL,HT(II))
     	      ENDIF
	   END IF
	END DO
C
C	Find and draw bonds.
C	If PWIDE < 0, fill bonds.
C	If PWIDE = 0, single line.
C	If PWIDE > 0, open bonds.
C
C	CALL NEWPEN(PWIDE(II))
	DWID  = 0.01
	NWID  = IFIX( 0.5 + 1.414**IABS(PWIDE(II)) )
	NWID1 = 1
	IF( PWIDE(II).EQ.0 ) NWID = 0
	IF( PWIDE(II).GT.0 ) THEN
	   DWID  = 0.01 * NWID
	   NWID  = 1
	   NWID1 = 2
	END IF
	DO J=1,JN-1
	   J1=J+1
	   DO J2=J1,JN
	      D=0.0
	      DO K=1,3
	         D=D+(XI(K,J)-XI(K,J2))**2
	         IF(D.GT.DMX2)GO TO 60
	      END DO
	      D=SQRT(D)
	      KJ1=KI(J)
	      KJ2=KI(J2)
	      IF(D.GT.RDEL(II)+RAD(II,KJ1)+RAD(II,KJ2))GO TO 60
	      IF(XI(3,J).LT.DZ3N.AND.XI(3,J2).LT.DZ3N)GO TO 60
	      IF(XI(3,J).GT.DZ3.AND.XI(3,J2).GT.DZ3)GO TO 60
	      XIN=XI(1,J)+XI(3,J)*STROX
	      YIN=XI(2,J)+XI(3,J)*STROY
	      XTR=XI(1,J2)+XI(3,J2)*STROX
	      YTR=XI(2,J2)+XI(3,J2)*STROY
	      IF((XIN.LT.XMPMIN).AND.(XTR.LT.XMPMIN))GO TO 60
	      IF((XIN.GT.XMPMAX).AND.(XTR.GT.XMPMAX))GO TO 60
	      IF((YIN.LT.YMPMIN).AND.(YTR.LT.YMPMIN))GO TO 60
	      IF((YIN.GT.YMPMAX).AND.(YTR.GT.YMPMAX))GO TO 60
	      IF( XI(3,J) .LT.DZ3N .OR. XI(3,J) .GT.DZ3 .OR.
	1         XI(3,J2).LT.DZ3N .OR. XI(3,J2).GT.DZ3 ) THEN
	         A1 = (XTR-XIN) / (XI(3,J2)-XI(3,J))
	         A2 = (YTR-YIN) / (XI(3,J2)-XI(3,J))
	         B1 = XIN - A1*XI(3,J)
	         B2 = YIN - A2*XI(3,J)
	         IF( XI(3,J).LT.DZ3N ) THEN
	            XIN = A1 * DZ3N + B1
	            YIN = A2 * DZ3N + B2
	         ELSE IF( XI(3,J).GT.DZ3 ) THEN
	            XIN = A1 * DZ3 + B1
	            YIN = A2 * DZ3 + B2
	         END IF
	         IF( XI(3,J2).LT.DZ3N ) THEN
	            XTR = A1 * DZ3N + B1
	            YTR = A2 * DZ3N + B2
	         ELSE IF( XI(3,J2).GT.DZ3 ) THEN
	            XTR = A1 * DZ3 + B1
	            YTR = A2 * DZ3 + B2
	         END IF
	      END IF
	      DDX = XTR - XIN 
	      DDY = YTR - YIN 
	      DDR = SQRT( DDX**2 + DDY**2 )
c	      write( 85, 111) j, j2, xin, yin, xtr, ytr
c111	      format( 2i5, 4f10.4)
	      IF( DDR.EQ.0.0 )THEN
	         DO ID = -NWID, NWID, NWID1
	            DDX  = DWID * ID
	            XIN1 = XIN - DDY
	            XTR1 = XTR - DDY
	            YIN1 = YIN - DDX
	            YTR1 = YTR - DDX
	            XIN1 = XIN - DDX
	            XTR1 = XTR - DDX
	            XIN2 = XIN + DDX
	            XTR2 = XTR + DDX
	            YIN1 = YIN - DDX
	            YTR1 = YTR + DDX
	            YIN2 = YIN - DDX
	            YTR2 = YTR + DDX
     		    CALL STRCON_LINEV1( XIN1, YIN1, XTR1, YTR1,
     +              XMPMIN, XMPMAX, YMPMIN, YMPMAX )
	            CALL STRCON_LINEV1( XIN2, YIN2, XTR2, YTR2,
     +              XMPMIN, XMPMAX, YMPMIN, YMPMAX )
	            CALL STRCON_LINEV1( XIN1, YIN1, XIN2, YIN2,
     +              XMPMIN, XMPMAX, YMPMIN, YMPMAX )
     		    CALL STRCON_LINEV1( XTR1, YTR1, XTR2, YTR2,
     +              XMPMIN, XMPMAX, YMPMIN, YMPMAX )
	         END DO
	      ELSE
	         DDX1  =  DWID * DDX / DDR
	         DDY1  =  DWID * DDY / DDR
	         DO ID = -NWID, NWID, NWID1
	            DDY  = DDY1 * ID
	            DDX  = DDX1 * ID
	            XIN1 = XIN - DDY
	            XTR1 = XTR - DDY
	            YIN1 = YIN + DDX
	            YTR1 = YTR + DDX
	            CALL STRCON_LINEV1( XIN1, YIN1, XTR1, YTR1,
     +              XMPMIN, XMPMAX, YMPMIN, YMPMAX )
	         END DO
	      END IF
 60	      CONTINUE
	   END DO
	END DO
C
 100	CONTINUE
	RETURN
	END

	SUBROUTINE STRCON_LINEV1( XA, YA, XB, YB, XMN, XMX, YMN, YMX)
	X1 = XA
	Y1 = YA
	X2 = XB
	Y2 = YB
	IF( X1.GE.XMN .AND. X1.LE.XMX .AND. Y1.GE.YMN .AND. Y1.LE.YMX .AND.
	1   X2.GE.XMN .AND. X2.LE.XMX .AND. Y2.GE.YMN .AND. Y2.LE.YMX ) THEN
	   CALL STRCON_LINEV( X1, Y1, X2, Y2 )
c	   iq = 0
c	   write( 85, 111) iq, x1, y1, x2, y2, xmn, xmx, ymn, ymx
c111	   format( i2, 8f9.4)
	   RETURN
	ELSE IF( ( X1.LE.XMN .AND. X2.LE.XMN ) .OR. 
	1        ( X1.GE.XMX .AND. X2.GE.XMX ) .OR. 
	2        ( Y1.LE.YMN .AND. Y2.LE.YMN ) .OR. 	
	3        ( Y1.GE.YMX .AND. Y2.GE.YMX ) )THEN
c	   iq = 2
c	   write( 85, 111) iq, x1, y1, x2, y2, xmn, xmx, ymn, ymx
	   RETURN
	ELSE
	   A1 = 0.0
	   IF( X2-X1.NE.0.0 ) A1 = ( Y2 - Y1 ) / ( X2 - X1 )
	   B1 = Y1 - A1 * X1 
	   IF( X1.LT.XMN ) THEN
	      X1 = XMN
	      Y1 = B1 + A1 * X1
	   END IF
	   IF( X2.LT.XMN ) THEN
	      X2 = XMN
	      Y2 = B1 + A1 * X2
	   END IF
	   IF( X1.GT.XMX ) THEN
	      X1 = XMX
	      Y1 = B1 + A1 * X1
	   END IF
	   IF( X2.GT.XMX ) THEN
	      X2 = XMX
	      Y2 = B1 + A1 * X2
	   END IF
	   A2 = 0.0
	   IF( Y2-Y1.NE.0.0 ) A2 = ( X2 - X1 ) / ( Y2 - Y1 )
	   B2 = X1 - A2 * Y1 
	   IF( Y1.LT.YMN ) THEN
	      Y1 = YMN
	      X1 = B2 + A2 * Y1
	   END IF
	   IF( Y2.LT.YMN ) THEN
	      Y2 = YMN
	      X2 = B2 + A2 * Y2
	   END IF
	   IF( Y1.GT.YMX ) THEN
	      Y1 = YMX
	      X1 = B2 + A2 * Y1
	   END IF
	   IF( Y2.GT.YMX ) THEN
	      Y2 = YMX
	      X2 = B2 + A2 * Y2
	   END IF
	   IF( X1.GE.XMN .AND. X1.LE.XMX .AND. Y1.GE.YMN .AND. Y1.LE.YMX .AND.
	1      X2.GE.XMN .AND. X2.LE.XMX .AND. Y2.GE.YMN .AND. Y2.LE.YMX ) THEN
c	      iq = 3
c	      write( 85, 111) iq, x1, y1, x2, y2, xmn, xmx, ymn, ymx
	      CALL STRCON_LINEV( X1, Y1, X2, Y2 )
	      RETURN
	   ELSE 
c	      iq = 4
c	      write( 85, 111) iq, x1, y1, x2, y2, xmn, xmx, ymn, ymx
	      RETURN
	   END IF
	END IF
	RETURN
	END

	SUBROUTINE STRCON_MODEL1(NDIR)
	COMMON /DAT/ TI(10), XMIN(3), XMAX(3), DX(3), XCHAR(3), ACHAR(3),
	1  IDIM(2),JDIM(2),NEXT,XSHX,XSHY,XSHZ,YSHY,YSHZ,ZSHZ,XMPMIN,
	2  XMPMAX, YMPMIN, YMPMAX, ISCALE, NSDONE, NIP, ZMAP, ABC,
	3  RMIN, RMAX, RINC, NSEC, ANG, RC
	COMMON /IOSET/ IO1,IO2,IO3,IO4,IO5,IO6,IO7,IO8
	COMMON/STRO/ISTRO,STROX,STROY,SSTROX,SSTROY,MSEC,MSEC2,IINZ
	COMMON /MOD/ AM(3,3,36),AT(3,36),XMIN3(3),XMAX3(3),
	1 DZ3,DZ3N,DZ2,DZ2N,DMX(4),CHAR(4,10),NSYM,PWIDE(4),ICHR(4),
	2 HT(4),HT2(4),SCALE,RAD(4,10),RDEL(4),IMASK(36,4),
	3 NSET(4,20),ATMFIL(36,4),VRAD(10),VTYP(10),IFVDW
	DIMENSION A(3,3),T(3),B(3,3),ABC(3)
	INTEGER*4 PWIDE,CHAR,VTYP
	LOGICAL*1 ATMFIL
	K1=1
	K2=2
	K3=3
	CALL STRCON_IROT(K1,K2,K3,NDIR)
	READ(IO2,1)NSYM
 1	FORMAT(I5)
	IF(NSYM.EQ.0)RETURN
	READ(IO2,13)IFVDW,(VTYP(I),VRAD(I),I=1,10)
 13	FORMAT(I5,10(I3,F5.1))
	DO 12 I=1,10
	IF(VTYP(I).LE.0.OR.VTYP(I).GT.4)VTYP(I)=0
	IF(VRAD(I).LT.0.0.OR.VRAD(I).GT.5.0)VRAD(I)=0.0
 12	CONTINUE
	DO 2 I=1,4
	DMX(I)=0.0
	PWIDE(I)=0
	ICHR(I)=0
	HT(I)=0.0
	READ(IO2,3)IQ,(ATMFIL(K,I),K=1,IQ)
	IF(IQ.LE.3)GO TO 2
 3	FORMAT(Q,32A1)
 	READ(IO2,11)DMX(I),PWIDE(I),ICHR(I),HT(I),(NSET(I,K),K=1,20)
 11	FORMAT(F5.2,2I5,F5.2,20I5)
	IF(NSET(I,2).EQ.0)NSET(I,2)=32000
	WRITE(IO6,111)I,(ATMFIL(K,I),K=1,IQ)
	WRITE(IO6,113)DMX(I),PWIDE(I),ICHR(I),
	1 HT(I),(NSET(I,K),K=1,20)
 111	FORMAT('1*****   INPUT PARAMETERS FOR MODEL PORTION',I2,/,
	1 '  MODEL FROM FILE ',32A1)
 113	FORMAT
	2('   INITIAL MAX BOND DIST.=',F6.2,/,'   PEN WIDTH =',I3,/,
	3 '   CHARACTER TYPE =',I3,/,'   CHARACTER HEIGHT =',F4.2,'MM.',
	4 10(/,'    ACCEPT ATOMS FROM #',I5,' TO #',I5))
	HT2(I) = HT(I)/2.0
	IF(ICHR(I).EQ.1)GO TO 32
	HT2(I)=0.0
	READ(IO2,33)(CHAR(I,J),J=1,10)
 33	FORMAT(10I4)
	IF(ICHR(I).NE.0)ICHR(I)=-1
	GO TO 34
 32	READ(IO2,35)(CHAR(I,J),J=1,10)
 35	FORMAT(10(3X,A1))
 34	CONTINUE
	READ(IO2,37)RDEL(I),(RAD(I,J),J=1,10)
 37	FORMAT(11F5.2)
	IF(ICHR(I))42,46,46
 42	WRITE(IO6,43)RDEL(I),(J,CHAR(I,J),RAD(I,J),VRAD(J),VTYP(J),
	1 J=1,10)
 43	FORMAT('0 INCREMENT IN BOND SEARCH =',F5.3,/,
	1'0 TYPE  CHARACTER  BOND RADIUS VDW RADIUS TYPE',/,
	2 (' ',I5,2X,I4,8X,F6.3,F10.3,I3))
	GO TO 48
 46	WRITE(IO6,47)RDEL(I),(J,CHAR(I,J),RAD(I,J),VRAD(J),VTYP(J),
	1 J=1,10)
 47	FORMAT('0 INCREMENT IN BOND SEARCH =',F5.3,/,
	1'0 TYPE  CHARACTER  BOND RADIUS VDW RADIUS TYPE',/,
	2 (' ',I5,5X,A1,8X,F6.3,F10.3,I3))
 48	CONTINUE
	READ(IO2,39)(IMASK(J,I),J=1,NSYM)
 39	FORMAT(36I2)
	WRITE(IO6,49)(IMASK(J,I),J=1,NSYM)
 49	FORMAT('0  SYMMETRY OPERATORS TO BE USED',36I2)
 2	CONTINUE
	DMXMX = 0.0
	DO K = 1, 4
	   IF( DMX(K).GT.DMXMX ) DMXMX = DMX(K)
	END DO
	IF( DMXMX.LT.3.5 ) DMXMX = 3.5
	XMIN3(1)=XMPMIN-DMXMX
	XMAX3(1)=XMPMAX+DMXMX
	XMIN3(2)=YMPMIN-DMXMX
	XMAX3(2)=YMPMAX+DMXMX
	DZ3=DX(3)*0.6666667*FLOAT(MSEC)
	DZ3N=-DZ3
	DZ2=DX(3)*0.5*FLOAT(MSEC)
	DZ2N=-DZ2
	DO 4 I=1,3
	DO 6 J=1,3
 6	B(I,J)=0.0
 4	CONTINUE
	B(1,1)=ABC(1)*XSHX
	B(2,2)=ABC(2)*YSHY
	B(3,3)=ABC(3)*ZSHZ
	B(1,2)=ABC(2)*XSHY
	B(1,3)=ABC(3)*XSHZ
	B(2,3)=ABC(3)*YSHZ
	IF(NSYM.LE.36)GO TO 22
	WRITE(IO6,213)NSYM
 213	FORMAT(' **** LIMIT OF 36 SYMMETRY OPERATORS EXCEDED, NSYM=',I3)
	STOP
 22	WRITE(IO6,21) NSYM
 21	FORMAT('1',I5,' SYMMETRY EQUIVALENT ATOMS TO BE GENERATED')
 	DO 10 I=1,NSYM
	READ(IO2,17)A(K1,1),A(K1,2),A(K1,3),A(K2,1),A(K2,2),
	1 A(K2,3),A(K3,1),A(K3,2),A(K3,3),T(K1),T(K2),T(K3)
	WRITE(IO6,23)(A(K,1),A(K,2),A(K,3),T(K),K=1,3)
 23	FORMAT('0',3(/,'  (',3F10.5,' ) +',F10.5))
 17	FORMAT(12F6.3)
	DO 14 I1=1,3
	AT(I1,I)=0.0
	DO 212 I2=1,3
	AM(I1,I2,I)=0.0
	DO 212 K=1,3
	AM(I1,I2,I)=AM(I1,I2,I)+B(I1,K)*A(K,I2)
 212	CONTINUE
	DO 14 K=1,3
	AT(I1,I)=AT(I1,I)+B(I1,K)*T(K)
 14	CONTINUE
 10	CONTINUE
	WRITE(IO6,27)
 27	FORMAT('1')
	RETURN
	END

     	SUBROUTINE STRCON_QPLOT(X,Y,IC)
C
C  PLOT A POINT OR A LINE OR MOVE PEN
C
C	X,Y=  COORDINATES, IN INCHES
C	IB= PEN CONTROL. IC=2, PEN IS DOWN DURING MOVE,
C			 IC=3, PEN IS UP DURING MOVE
C			 IC=-2 OR -3, A NEW ORIGIN IS DEFINED
C
     	PARAMETER CIN= 25.4
	INTEGER   IC, IB, IP
	COMMON/OFFSET/XOFF,XFAC,YOFF,YFAC

	XN = X
	YN = Y
	IB = IC
     	IP = IABS(IB)
     	IF( IP.GT.10 ) THEN
     	    IB = ISIGN(IP-10,IB)
     	    XN = (XN-XOFF) / XFAC / CIN
     	    YN = (YN-YOFF) / YFAC / CIN
     	ENDIF
	CALL IFPLOT_VECT( XN, YN, IB )
	RETURN
	END

	SUBROUTINE STRCON_OFFSET(XOFFX,XFACX,YOFFX,YFACX)
C
	COMMON/OFFSET/XOFF,XFAC,YOFF,YFAC
	DATA XOFF,XFAC,YOFF,YFAC /0.0,1.0,0.0,1.0/
C
	XOFF = XOFFX
	XFAC = XFACX
	YOFF = YOFFX
	YFAC = YFACX
	RETURN
	END

	SUBROUTINE STRCON_DOTLNV(X1,Y1,X2,Y2,DOT,SP)
C
C	SUBROUTINE TO DRAW A DOTTED LINE FROM X1,Y1 TO X2,Y2
C	WITH DOTS OF LENGTH 'DOT' AND SEPERATED BY SPACINGS 'SP'.
C
C
	XI=X1
	YI=Y1
	X0=XI
	Y0=YI
	XT=X2
	YT=Y2
	CALL STRCON_QPLOT(XI,YI,+13)
	DX=XT-XI
	DY=YT-YI
	D=SQRT(DX**2+DY**2)
	N=D/(DOT+SP)
	IF(N.EQ.0)GO TO 10
	DXS=SP*DX/D
	DYS=SP*DY/D
	DXD=DOT*DX/D
	DYD=DOT*DY/D
	X0=XI+DXD/2.0
	Y0=YI+DYD/2.0
	DO 20 I=1,N
	CALL STRCON_QPLOT(X0,Y0,+12)
	X3=X0+DXS
	Y3=Y0+DYS
	CALL STRCON_QPLOT(X3,Y3,+13)
	X0=X3+DXD
	Y0=Y3+DYD
 20	CONTINUE
	X0=X3+DXD/2.0
	Y0=Y3+DYD/2.0
	CALL STRCON_QPLOT(X0,Y0,+12)
	DX=XT-X0
	DY=YT-Y0
	D=SQRT(DX**2+DY**2)
 10	IF(D.EQ.0)RETURN
	DX=DX/6.0
	DY=DY/6.0
 	CALL STRCON_QPLOT(X0+DX,Y0+DY,+12)
	CALL STRCON_QPLOT(XT-DX,YT-DY,+13)
	CALL STRCON_QPLOT(XT,YT,+12)
	RETURN
	END

	SUBROUTINE STRCON_TACON(JDIM1,JDIM2,RMIN,RMAX,RINC)
C
C      tacon, triangular approximation contouring, approximates the
C      surface to be contoured as a connected set of triangular planes
C      and derives the contours in piecewise parametric form from the
C      triangular planes.
C
	COMMON/MAPDAT/  XMAP(150), YMAP(150), XMAPY(150)
	COMMON /IOSET/ IO1,IO2,IO3,IO4,IO5,IO6,IO7,IO8
	DIMENSION RV(5), XV(5), YV(5), RT(4), XT(4), YT(4), IADD(2)
	LOGICAL IVINT(3)
	REAL*4 IX(3),IY(3)
	INTEGER*2 RMAP(2,150),RTST,IMAP(2,152)
	EQUIVALENCE (RMAP(1,1),IMAP(1,2))
	DATA RTST/32767/
	IMAX = JDIM1 - 1
	JMAX = JDIM2 - 1
	DX=XMAP(2)-XMAP(1)
	REWIND IO3
	READ(IO3)(RMAP(2,I),I=1,JDIM1)
C      START THROUGH THE GRID
	DO 51 J = 1, JMAX
	XMY=XMAPY(J)
	XMY1=XMAPY(J+1)
	ID=0
	IF(ABS(XMY1-XMY)/DX.LT.0.5)GO TO 152
	IF(XMY.GT.XMY1)GO TO 154
	IMAP(2,JDIM1+2)=RTST
	ID=1
	XMY=XMY+DX
	GO TO 152
 154	IMAP(2,1)=RTST
	ID=-1
	XMY=XMY-DX
 152	CONTINUE
	DO 151 I=1,JDIM1
	   I1=I+ID
 151	   RMAP(1,I)=RMAP(2,I1)
	YV(1) = YMAP(J)
	YV(2) = YMAP(J)
	YV(3) = YMAP(J+1)
	YV(4) = YMAP(J+1)
	READ(IO3)(RMAP(2,I),I=1,JDIM1)
	DO 52 I = 1, IMAX
C      GET A RECTANGLE WITH CORNERS ON THE GRID POINTS
	IF((RMAP(2,I).EQ.RTST).OR.(RMAP(1,I).EQ.RTST).OR.
	1  (RMAP(2,I+1).EQ.RTST).OR.(RMAP(1,I+1).EQ.RTST))GO TO 52
	RV(1) = RMAP(1, I)
	XV(1) = XMAP(I) + XMY
	RV(2) = RMAP(1, I+1)
	XV(2) = XMAP(I+1) + XMY
	RV(3) = RMAP(2, I+1)
	XV(3) = XMAP(I+1) + XMY1
	RV(4) = RMAP(2, I)
	XV(4) = XMAP(I) + XMY1
C      FIND OUT HOW MANY, IF ANY, CONTOURS ARE CROSSED WITHIN THE REC-
C      TANGLE.  SKIP OR SET UP LIMITS ACCORDINGLY
	RMPMAX = AMAX1(RV(1), RV(2), RV(3), RV(4))
	IF(RMPMAX .LT. RMIN) GO TO 52
	RMPMIN = AMIN1(RV(1), RV(2), RV(3), RV(4))
	IF(RMPMIN .GT. RMAX) GO TO 52
	MSTART = (RMPMIN - RMIN)/RINC + 1.01
	IF(MSTART .LE. 0) MSTART = 1
	IF(RMPMAX .GT. RMAX) RMPMAX = RMAX
	MSTOP = (RMPMAX - RMIN)/RINC + 1.01
	IF(MSTOP .GE. MSTART) GO TO 60
	WRITE(IO6, 101) MSTART, MSTOP, I, J, RMPMAX, RMPMIN, (RV(II),
	1  II = 1, 4)
 101	FORMAT('0TACON ERROR'/' ', 4(I15, 5X)/' ', 1P6E10.5)
	CALL EXIT
 60	CONTINUE
	RV(5) = RV(1)
	XV(5) = XV(1)
	YV(5) = YV(1)
C      GET THE COMMON VERTEX OF ALL FOUR TRIANGLES WITHIN THIS RECTANGLE
C      BY INTERPOLATING A VALUE AT THE CENTER OF THE RECTANGLE
	RT(1) = (RV(1) + RV(2) + RV(3) + RV(4))/4.
	XT(1) = (XV(1) + XV(3))/2.
	YT(1) = (YV(1) + YV(3))/2.
C      START AROUND THE RECTANGLE
	DO 53 K = 1, 4
C      GET A TRIANGLE
	RT(2) = RV(K)
	XT(2) = XV(K)
	YT(2) = YV(K)
	RT(3) = RV(K+1)
	XT(3) = XV(K+1)
	YT(3) = YV(K+1)
	RT(4) = RT(1)
	XT(4) = XT(1)
	YT(4) = YT(1)
C      CONTOUR THE TRIANGLE
	DO 54 M = MSTART, MSTOP
	XM = M - 1
	RLVL = RMIN + RINC*XM
	DO 9 L = 1, 3
 9	IVINT(L) = .FALSE.
	NVINT = 1
	DO 10 L = 1, 3
	IF(RT(L) .NE. RLVL) GO TO 10
	NVINT = NVINT + 1
	IVINT(L) = .TRUE.
 10	CONTINUE
	GO TO (20, 30, 40, 54), NVINT
C      CONTOUR LEVEL DOES NOT INTERSECT ANY VERTICES
 20	CONTINUE
	NCROSS = 0
	DO 25 L = 1, 3
	IF(RT(L) .GT. RT(L+1)) GO TO 21
C      RT(L + IADD(2)) IS THE LARGER R VALUE
	IADD(1) = 0
	IADD(2) = 1
	GO TO 22
 21	CONTINUE
	IADD(1) = 1
	IADD(2) = 0
 22	CONTINUE
	IF((RLVL .LE. RT(L+IADD(1)) ).OR.(RLVL .GE. RT(L+IADD(2)) ))
	1  GO TO 25
	NCROSS = NCROSS + 1
	DELTA = (RLVL - RT(L+IADD(1)))/(RT(L+IADD(2)) - RT(L+IADD(1)))
	X = (XT(L+IADD(2))-XT(L+IADD(1)))*DELTA + XT(L+IADD(1))
	IX(NCROSS) = X
	Y = (YT(L+IADD(2)) - YT(L+IADD(1)))*DELTA + YT(L+IADD(1))
	IY(NCROSS) = Y
 25	CONTINUE
	IF(NCROSS .EQ. 0) GO TO 54
	IF(NCROSS .EQ. 2) GO TO 26
	WRITE(IO6, 120) NCROSS, I, J, RLVL, ((RT(II),XT(II),YT(II)),II=
	1  1, 3)
 120	FORMAT('0ERROR>>>>>>>>>>NCROSS = ', I1/' I = ', I3, ', J = ',
	1  I3, ', RLVL = ', F9.4/'0RT', 8X, 'XT', 8X, 'YT'/
	2  (' ', 3(F9.4, 1X), 1X))
	GO TO 54
26  	CONTINUE
	CALL STRCON_LINEV(IX(1), IY(1), IX(2), IY(2))
	GO TO 54
C      CONTOUR LEVEL INTERSECTS ONE VERTEX
 30	CONTINUE
	DO 31 L = 1, 3
 31	IF(IVINT(L)) IV = L
	GO TO (32, 33, 34), IV
 32	IV1 = 2
	IV2 = 3
	GO TO 35
 33	IV1 = 1
	IV2 = 3
	GO TO 35
 34	IV1 = 1
	IV2 = 2
 35	CONTINUE
	IF(RT(IV2) .GT. RT(IV1)) GO TO 36
	IV3 = IV2
	IV2 = IV1
	IV1 = IV3
 36	CONTINUE
	IF((RLVL .LE. RT(IV1)).OR.(RLVL.GE.RT(IV2))) GO TO 54
	DELTA = (RLVL - RT(IV1))/(RT(IV2) - RT(IV1))
	X = (XT(IV2) - XT(IV1))*DELTA + XT(IV1)
	IX(1) = X
	Y = (YT(IV2) - YT(IV1))*DELTA + YT(IV1)
	IY(1) = Y
	IX(2) = XT(IV)
	IY(2) = YT(IV)
	CALL STRCON_LINEV(IX(1), IY(1), IX(2), IY(2))
	GO TO 54
C      CONTOUR LEVEL INTERSECTS TWO VERTICES
 40	CONTINUE
	II = 1
	DO 41 L = 1, 3
	IF(.NOT.IVINT(L)) GO TO 41
	IX(II) = XT(L)
	IY(II) = YT(L)
	II = II + 1
 41	CONTINUE
	CALL STRCON_LINEV(IX(1), IY(1), IX(2), IY(2))
 54	CONTINUE
 53	CONTINUE
 52	CONTINUE
 51	CONTINUE
	RETURN
	END

	SUBROUTINE STRCON_TACOF(JDIM1,JDIM2,RMIN,RMAX,RINC)
C      TACON, TRIANGULAR APPROXIMATION CONTOURING, APPROXIMATES THE
C      SURFACE TO BE CONTOURED AS A CONNECTED SET OF TRIANGULAR PLANES
C      AND DERIVES THE CONTOURS IN PIECEWISE PARAMETRIC FORM FROM THE
C      TRIANGULAR PLANES.
	COMMON/MAPDAT/  XMAP(150), YMAP(150), XMAPY(150)
	COMMON /IOSET/ IO1,IO2,IO3,IO4,IO5,IO6,IO7,IO8
	DIMENSION RV(5), XV(5), YV(5), RT(4), XT(4), YT(4), IADD(2)
	LOGICAL IVINT(3)
	REAL*4 IX(3),IY(3)
	INTEGER*2 RMAP(2,150),RTST,IMAP(2,152)
	EQUIVALENCE (RMAP(1,1),IMAP(1,2))
	DATA RTST/32767/
	DOT=0.15
	SP=0.3
	IMAX = JDIM1 - 1
	JMAX = JDIM2 - 1
	DX=XMAP(2)-XMAP(1)
	REWIND IO3
	READ(IO3)(RMAP(2,I),I=1,JDIM1)
C      START THROUGH THE GRID
	DO 51 J = 1, JMAX
	XMY=XMAPY(J)
	XMY1=XMAPY(J+1)
	ID=0
	IF(XMY*XMY1.GE.0.0)GO TO 152
	IF(ABS(XMY1-XMY)/DX.LT.0.5)GO TO 152
	IF(XMY.GT.XMY1)GO TO 154
	IMAP(2,JDIM1+2)=RTST
	ID=1
	XMY=XMY+DX
	GO TO 152
 154	IMAP(2,1)=RTST
	ID=-1
	XMY=XMY-DX
 152	CONTINUE
	DO 151 I=1,JDIM1
	   I1=I+ID
 151	   RMAP(1,I)=RMAP(2,I1)
	YV(1) = YMAP(J)
	YV(2) = YMAP(J)
	YV(3) = YMAP(J+1)
	YV(4) = YMAP(J+1)
	READ(IO3)(RMAP(2,I),I=1,JDIM1)
	DO 52 I = 1, IMAX
	IF((RMAP(2,I).EQ.RTST).OR.(RMAP(1,I).EQ.RTST).OR.
	1  (RMAP(2,I+1).EQ.RTST).OR.(RMAP(1,I+1).EQ.RTST))GO TO 52
C      GET A RECTANGLE WITH CORNERS ON THE GRID POINTS
	RV(1) = -RMAP(1, I)
	XV(1) = XMAP(I) + XMY
	RV(2) = -RMAP(1, I+1)
	XV(2) = XMAP(I+1) + XMY
	RV(3) = -RMAP(2, I+1)
	XV(3) = XMAP(I+1) + XMY1
	RV(4) = -RMAP(2, I)
	XV(4) = XMAP(I) + XMY1
C      FIND OUT HOW MANY, IF ANY, CONTOURS ARE CROSSED WITHIN THE REC-
C      TANGLE.  SKIP OR SET UP LIMITS ACCORDINGLY
	RMPMAX = AMAX1(RV(1), RV(2), RV(3), RV(4))
	IF(RMPMAX .LT. RMIN) GO TO 52
	RMPMIN = AMIN1(RV(1), RV(2), RV(3), RV(4))
	IF(RMPMIN .GT. RMAX) GO TO 52
	MSTART = (RMPMIN - RMIN)/RINC + 1.01
	IF(MSTART .LE. 0) MSTART = 1
	IF(RMPMAX .GT. RMAX) RMPMAX = RMAX
	MSTOP = (RMPMAX - RMIN)/RINC + 1.01
	IF(MSTOP .GE. MSTART) GO TO 60
	WRITE(IO6, 101) MSTART, MSTOP, I, J, RMPMAX, RMPMIN, (RV(II),
	1  II = 1, 4)
 101	FORMAT('0TACON ERROR'/' ', 4(I15, 5X)/' ', 1P6E10.5)
	CALL EXIT
 60	CONTINUE
	RV(5) = RV(1)
	XV(5) = XV(1)
	YV(5) = YV(1)
C      GET THE COMMON VERTEX OF ALL FOUR TRIANGLES WITHIN THIS RECTANGLE
C      BY INTERPOLATING A VALUE AT THE CENTER OF THE RECTANGLE
	RT(1) = (RV(1) + RV(2) + RV(3) + RV(4))/4.
	XT(1) = (XV(1) + XV(3))/2.
	YT(1) = (YV(1) + YV(3))/2.
C      START AROUND THE RECTANGLE
	DO 53 K = 1, 4
C      GET A TRIANGLE
	RT(2) = RV(K)
	XT(2) = XV(K)
	YT(2) = YV(K)
	RT(3) = RV(K+1)
	XT(3) = XV(K+1)
	YT(3) = YV(K+1)
	RT(4) = RT(1)
	XT(4) = XT(1)
	YT(4) = YT(1)
C      CONTOUR THE TRIANGLE
	DO 54 M = MSTART, MSTOP
	XM = M - 1
	RLVL = RMIN + RINC*XM
	DO 9 L = 1, 3
 9	IVINT(L) = .FALSE.
	NVINT = 1
	DO 10 L = 1, 3
	IF(RT(L) .NE. RLVL) GO TO 10
	NVINT = NVINT + 1
	IVINT(L) = .TRUE.
 10	CONTINUE
	GO TO (20, 30, 40, 54), NVINT
C      CONTOUR LEVEL DOES NOT INTERSECT ANY VERTICES
 20	CONTINUE
	NCROSS = 0
	DO 25 L = 1, 3
	IF(RT(L) .GT. RT(L+1)) GO TO 21
C      RT(L + IADD(2)) IS THE LARGER R VALUE
	IADD(1) = 0
	IADD(2) = 1
	GO TO 22
 21	CONTINUE
	IADD(1) = 1
	IADD(2) = 0
 22	CONTINUE
	IF((RLVL .LE. RT(L+IADD(1)) ).OR.(RLVL .GE. RT(L+IADD(2)) ))
	1  GO TO 25
	NCROSS = NCROSS + 1
	DELTA = (RLVL - RT(L+IADD(1)))/(RT(L+IADD(2)) - RT(L+IADD(1)))
	X = (XT(L+IADD(2))-XT(L+IADD(1)))*DELTA + XT(L+IADD(1))
	IX(NCROSS) = X
	Y = (YT(L+IADD(2)) - YT(L+IADD(1)))*DELTA + YT(L+IADD(1))
	IY(NCROSS) = Y
 25	CONTINUE
	IF(NCROSS .EQ. 0) GO TO 54
	IF(NCROSS .EQ. 2) GO TO 26
	WRITE(IO6, 120) NCROSS, I, J, RLVL, ((RT(II),XT(II),YT(II)),II=
	1  1, 3)
 120	FORMAT('0ERROR>>>>>>>>>>NCROSS = ', I1/' I = ', I3, ', J = ',
	1  I3, ', RLVL = ', F9.4/'0RT', 8X, 'XT', 8X, 'YT'/
	2  (' ', 3(F9.4, 1X), 1X))
	GO TO 54
26  	CONTINUE
	CALL STRCON_DOTLNV(IX(1), IY(1), IX(2), IY(2),DOT,SP)
	GO TO 54
C      CONTOUR LEVEL INTERSECTS ONE VERTEX
 30	CONTINUE
	DO 31 L = 1, 3
 31	IF(IVINT(L)) IV = L
	GO TO (32, 33, 34), IV
 32	IV1 = 2
	IV2 = 3
	GO TO 35
 33	IV1 = 1
	IV2 = 3
	GO TO 35
 34	IV1 = 1
	IV2 = 2
 35	CONTINUE
	IF(RT(IV2) .GT. RT(IV1)) GO TO 36
	IV3 = IV2
	IV2 = IV1
	IV1 = IV3
 36	CONTINUE
	IF((RLVL .LE. RT(IV1)).OR.(RLVL.GE.RT(IV2))) GO TO 54
	DELTA = (RLVL - RT(IV1))/(RT(IV2) - RT(IV1))
	X = (XT(IV2) - XT(IV1))*DELTA + XT(IV1)
	IX(1) = X
	Y = (YT(IV2) - YT(IV1))*DELTA + YT(IV1)
	IY(1) = Y
	IX(2) = XT(IV)
	IY(2) = YT(IV)
	CALL STRCON_DOTLNV(IX(1), IY(1), IX(2), IY(2),DOT,SP)
	GO TO 54
C      CONTOUR LEVEL INTERSECTS TWO VERTICES
 40	CONTINUE
	II = 1
	DO 41 L = 1, 3
	IF(.NOT.IVINT(L)) GO TO 41
	IX(II) = XT(L)
	IY(II) = YT(L)
	II = II + 1
 41	CONTINUE
	CALL STRCON_DOTLNV(IX(1), IY(1), IX(2), IY(2),DOT,SP)
 54	CONTINUE
 53	CONTINUE
 52	CONTINUE
 51	CONTINUE
	RETURN
	END

	SUBROUTINE STRCON_HACON(JDIM1,JDIM2,RMIN,RMAX,RINC)
C      TACON, TRIANGULAR APPROXIMATION CONTOURING, APPROXIMATES THE
C      SURFACE TO BE CONTOURED AS A CONNECTED SET OF TRIANGULAR PLANES
C      AND DERIVES THE CONTOURS IN PIECEWISE PARAMETRIC FORM FROM THE
C      TRIANGULAR PLANES.
	COMMON/MAPDAT/  XMAP(150), YMAP(150), XMAPY(150)
	COMMON /IOSET/ IO1,IO2,IO3,IO4,IO5,IO6,IO7,IO8
	DIMENSION RV(5), XV(5), YV(5), RT(4), XT(4), YT(4), IADD(2)
	LOGICAL IVINT(3)
	REAL*4 IX(3),IY(3)
	INTEGER*2 RMAP(2,150),RTST,IMAP(2,152)
	EQUIVALENCE (RMAP(1,1),IMAP(1,2))
	DATA RTST/32767/
	IMAX = JDIM1 - 1
	JMAX = JDIM2 - 1
	DX=XMAP(2)-XMAP(1)
	REWIND IO3
	READ(IO3)(RMAP(2,I),I=1,JDIM1)
C      START THROUGH THE GRID
	DO 51 J = 1, JMAX
	XMY=XMAPY(J)
	XMY1=XMAPY(J+1)
	IF(XMY1.GT.XMY)GO TO 154
	J1 = 2
	J2 = 1
	J3 = 3
	J4 = 4
	GO TO 152
 154	J1 = 1
	J2 = 2
	J3 = 4
	J4 = 3
 152	CONTINUE
	DO 151 I=1,JDIM1
 151	   RMAP(1,I)=RMAP(2,I)
	YV(1) = YMAP(J)
	YV(2) = YMAP(J)
	YV(3) = YMAP(J+1)
	YV(4) = YMAP(J+1)
	READ(IO3)(RMAP(2,I),I=1,JDIM1)
	DO 52 I = 1, IMAX
C      GET A RECTANGLE WITH CORNERS ON THE GRID POINTS
	IF((RMAP(2,I).EQ.RTST).OR.(RMAP(1,I).EQ.RTST).OR.
	1  (RMAP(2,I+1).EQ.RTST).OR.(RMAP(1,I+1).EQ.RTST))GO TO 52
	RV(1) = RMAP(1, I)
	XV(1) = XMAP(I) + XMY
	RV(2) = RMAP(1, I+1)
	XV(2) = XMAP(I+1) + XMY
	RV(3) = RMAP(2, I+1)
	XV(3) = XMAP(I+1) + XMY1
	RV(4) = RMAP(2, I)
	XV(4) = XMAP(I) + XMY1
C      FIND OUT HOW MANY, IF ANY, CONTOURS ARE CROSSED WITHIN THE REC-
C      TANGLE.  SKIP OR SET UP LIMITS ACCORDINGLY
	RMPMAX = AMAX1(RV(1), RV(2), RV(3), RV(4))
	IF(RMPMAX .LT. RMIN) GO TO 52
	RMPMIN = AMIN1(RV(1), RV(2), RV(3), RV(4))
	IF(RMPMIN .GT. RMAX) GO TO 52
	MSTART = (RMPMIN - RMIN)/RINC + 1.01
	IF(MSTART .LE. 0) MSTART = 1
	IF(RMPMAX .GT. RMAX) RMPMAX = RMAX
	MSTOP = (RMPMAX - RMIN)/RINC + 1.01
	IF(MSTOP .GE. MSTART) GO TO 60
	WRITE(IO6, 101) MSTART, MSTOP, I, J, RMPMAX, RMPMIN, (RV(II),
	1  II = 1, 4)
 101	FORMAT('0TACON ERROR'/' ', 4(I15, 5X)/' ', 1P6E10.5)
	CALL EXIT
 60	CONTINUE
C      START AROUND THE RECTANGLE
	DO 53 K = 1, 2
C      GET A TRIANGLE
	IF(K.EQ.2)GO TO 8
	RT(1)=RV(J1)
	XT(1)=XV(J1)
	YT(1)=YV(J1)
	RT(2)=RV(J2)
	XT(2)=XV(J2)
	YT(2)=YV(J2)
	RT(3)=RV(J3)
	XT(3)=XV(J3)
	YT(3)=YV(J3)
	GO TO 7
 8	RT(1)=RV(J2)
	XT(1)=XV(J2)
	YT(1)=YV(J2)
	RT(2)=RV(J3)
	XT(2)=XV(J3)
	YT(2)=YV(J3)
	RT(3)=RV(J4)
	XT(3)=XV(J4)
	YT(3)=YV(J4)
 7	RT(4)=RT(1)
	XT(4)=XT(1)
	YT(4)=YT(1)
C      CONTOUR THE TRIANGLE
	DO 54 M = MSTART, MSTOP
	XM = M - 1
	RLVL = RMIN + RINC*XM
	DO 9 L = 1, 3
 9	IVINT(L) = .FALSE.
	NVINT = 1
	DO 10 L = 1, 3
	IF(RT(L) .NE. RLVL) GO TO 10
	NVINT = NVINT + 1
	IVINT(L) = .TRUE.
 10	CONTINUE
	GO TO (20, 30, 40, 54), NVINT
C      CONTOUR LEVEL DOES NOT INTERSECT ANY VERTICES
 20	CONTINUE
	NCROSS = 0
	DO 25 L = 1, 3
	IF(RT(L) .GT. RT(L+1)) GO TO 21
C      RT(L + IADD(2)) IS THE LARGER R VALUE
	IADD(1) = 0
	IADD(2) = 1
	GO TO 22
 21	CONTINUE
	IADD(1) = 1
	IADD(2) = 0
 22	CONTINUE
	IF((RLVL .LE. RT(L+IADD(1)) ).OR.(RLVL .GE. RT(L+IADD(2)) ))
	1  GO TO 25
	NCROSS = NCROSS + 1
	DELTA = (RLVL - RT(L+IADD(1)))/(RT(L+IADD(2)) - RT(L+IADD(1)))
	X = (XT(L+IADD(2))-XT(L+IADD(1)))*DELTA + XT(L+IADD(1))
	IX(NCROSS) = X
	Y = (YT(L+IADD(2)) - YT(L+IADD(1)))*DELTA + YT(L+IADD(1))
	IY(NCROSS) = Y
 25	CONTINUE
	IF(NCROSS .EQ. 0) GO TO 54
	IF(NCROSS .EQ. 2) GO TO 26
	WRITE(IO6, 120) NCROSS, I, J, RLVL, ((RT(II),XT(II),YT(II)),II=
	1  1, 3)
 120	FORMAT('0ERROR>>>>>>>>>>NCROSS = ', I1/' I = ', I3, ', J = ',
	1  I3, ', RLVL = ', F9.4/'0RT', 8X, 'XT', 8X, 'YT'/
	1  (' ', 3(F9.4, 1X), 1X))
	GO TO 54
26  	CONTINUE
	CALL STRCON_LINEV(IX(1), IY(1), IX(2), IY(2))
	GO TO 54
C      CONTOUR LEVEL INTERSECTS ONE VERTEX
 30	CONTINUE
	DO 31 L = 1, 3
 31	IF(IVINT(L)) IV = L
	GO TO (32, 33, 34), IV
 32	IV1 = 2
	IV2 = 3
	GO TO 35
 33	IV1 = 1
	IV2 = 3
	GO TO 35
 34	IV1 = 1
	IV2 = 2
 35	CONTINUE
	IF(RT(IV2) .GT. RT(IV1)) GO TO 36
	IV3 = IV2
	IV2 = IV1
	IV1 = IV3
 36	CONTINUE
	IF((RLVL .LE. RT(IV1)).OR.(RLVL.GE.RT(IV2))) GO TO 54
	DELTA = (RLVL - RT(IV1))/(RT(IV2) - RT(IV1))
	X = (XT(IV2) - XT(IV1))*DELTA + XT(IV1)
	IX(1) = X
	Y = (YT(IV2) - YT(IV1))*DELTA + YT(IV1)
	IY(1) = Y
	IX(2) = XT(IV)
	IY(2) = YT(IV)
	CALL STRCON_LINEV(IX(1), IY(1), IX(2), IY(2))
	GO TO 54
C      CONTOUR LEVEL INTERSECTS TWO VERTICES
 40	CONTINUE
	II = 1
	DO 41 L = 1, 3
	IF(.NOT.IVINT(L)) GO TO 41
	IX(II) = XT(L)
	IY(II) = YT(L)
	II = II + 1
 41	CONTINUE
	CALL STRCON_LINEV(IX(1), IY(1), IX(2), IY(2))
 54	CONTINUE
 53	CONTINUE
 52	CONTINUE
 51	CONTINUE
	RETURN
	END

	SUBROUTINE STRCON_HACOF(JDIM1,JDIM2,RMIN,RMAX,RINC)
C      TACON, TRIANGULAR APPROXIMATION CONTOURING, APPROXIMATES THE
C      SURFACE TO BE CONTOURED AS A CONNECTED SET OF TRIANGULAR PLANES
C      AND DERIVES THE CONTOURS IN PIECEWISE PARAMETRIC FORM FROM THE
C      TRIANGULAR PLANES.
	COMMON/MAPDAT/  XMAP(150), YMAP(150), XMAPY(150)
	COMMON /IOSET/ IO1,IO2,IO3,IO4,IO5,IO6,IO7,IO8
	DIMENSION RV(5), XV(5), YV(5), RT(4), XT(4), YT(4), IADD(2)
	LOGICAL IVINT(3)
	REAL*4 IX(3),IY(3)
	INTEGER*2 RMAP(2,150),RTST,IMAP(2,152)
	EQUIVALENCE (RMAP(1,1),IMAP(1,2))
	DATA RTST/32767/
	DOT=0.15
	SP=0.3
	IMAX = JDIM1 - 1
	JMAX = JDIM2 - 1
	DX=XMAP(2)-XMAP(1)
	REWIND IO3
	READ(IO3)(RMAP(2,I),I=1,JDIM1)
C      START THROUGH THE GRID
	DO 51 J = 1, JMAX
	XMY=XMAPY(J)
	XMY1=XMAPY(J+1)
	IF(XMY1.GT.XMY)GO TO 154
	J1 = 2
	J2 = 1
	J3 = 3
	J4 = 4
	GO TO 152
 154	J1 = 1
	J2 = 2
	J3 = 4
	J4 = 3
 152	CONTINUE
	DO 151 I=1,JDIM1
 151	   RMAP(1,I)=RMAP(2,I)
	YV(1) = YMAP(J)
	YV(2) = YMAP(J)
	YV(3) = YMAP(J+1)
	YV(4) = YMAP(J+1)
	READ(IO3)(RMAP(2,I),I=1,JDIM1)
	DO 52 I = 1, IMAX
	IF((RMAP(2,I).EQ.RTST).OR.(RMAP(1,I).EQ.RTST).OR.
	1  (RMAP(2,I+1).EQ.RTST).OR.(RMAP(1,I+1).EQ.RTST))GO TO 52
C      GET A RECTANGLE WITH CORNERS ON THE GRID POINTS
	RV(1) = -RMAP(1, I)
	XV(1) = XMAP(I) + XMY
	RV(2) = -RMAP(1, I+1)
	XV(2) = XMAP(I+1) + XMY
	RV(3) = -RMAP(2, I+1)
	XV(3) = XMAP(I+1) + XMY1
	RV(4) = -RMAP(2, I)
	XV(4) = XMAP(I) + XMY1
C      FIND OUT HOW MANY, IF ANY, CONTOURS ARE CROSSED WITHIN THE REC-
C      TANGLE.  SKIP OR SET UP LIMITS ACCORDINGLY
	RMPMAX = AMAX1(RV(1), RV(2), RV(3), RV(4))
	IF(RMPMAX .LT. RMIN) GO TO 52
	RMPMIN = AMIN1(RV(1), RV(2), RV(3), RV(4))
	IF(RMPMIN .GT. RMAX) GO TO 52
	MSTART = (RMPMIN - RMIN)/RINC + 1.01
	IF(MSTART .LE. 0) MSTART = 1
	IF(RMPMAX .GT. RMAX) RMPMAX = RMAX
	MSTOP = (RMPMAX - RMIN)/RINC + 1.01
	IF(MSTOP .GE. MSTART) GO TO 60
	WRITE(IO6, 101) MSTART, MSTOP, I, J, RMPMAX, RMPMIN, (RV(II),
	1  II = 1, 4)
 101	FORMAT('0TACON ERROR'/' ', 4(I15, 5X)/' ', 1P6E10.5)
	CALL EXIT
 60	CONTINUE
C      START AROUND THE RECTANGLE
	DO 53 K = 1, 2
C      GET A TRIANGLE
	IF(K.EQ.2)GO TO 8
	RT(1)=RV(J1)
	XT(1)=XV(J1)
	YT(1)=YV(J1)
	RT(2)=RV(J2)
	XT(2)=XV(J2)
	YT(2)=YV(J2)
	RT(3)=RV(J3)
	XT(3)=XV(J3)
	YT(3)=YV(J3)
	GO TO 7
 8	RT(1)=RV(J2)
	XT(1)=XV(J2)
	YT(1)=YV(J2)
	RT(2)=RV(J3)
	XT(2)=XV(J3)
	YT(2)=YV(J3)
	RT(3)=RV(J4)
	XT(3)=XV(J4)
	YT(3)=YV(J4)
 7	RT(4)=RT(1)
	XT(4)=XT(1)
	YT(4)=YT(1)
C      CONTOUR THE TRIANGLE
	DO 54 M = MSTART, MSTOP
	XM = M - 1
	RLVL = RMIN + RINC*XM
	DO 9 L = 1, 3
 9	IVINT(L) = .FALSE.
	NVINT = 1
	DO 10 L = 1, 3
	IF(RT(L) .NE. RLVL) GO TO 10
	NVINT = NVINT + 1
	IVINT(L) = .TRUE.
 10	CONTINUE
	GO TO (20, 30, 40, 54), NVINT
C      CONTOUR LEVEL DOES NOT INTERSECT ANY VERTICES
 20	CONTINUE
	NCROSS = 0
	DO 25 L = 1, 3
	IF(RT(L) .GT. RT(L+1)) GO TO 21
C      RT(L + IADD(2)) IS THE LARGER R VALUE
	IADD(1) = 0
	IADD(2) = 1
	GO TO 22
 21	CONTINUE
	IADD(1) = 1
	IADD(2) = 0
 22	CONTINUE
	IF((RLVL .LE. RT(L+IADD(1)) ).OR.(RLVL .GE. RT(L+IADD(2)) ))
	1  GO TO 25
	NCROSS = NCROSS + 1
	DELTA = (RLVL - RT(L+IADD(1)))/(RT(L+IADD(2)) - RT(L+IADD(1)))
	X = (XT(L+IADD(2))-XT(L+IADD(1)))*DELTA + XT(L+IADD(1))
	IX(NCROSS) = X
	Y = (YT(L+IADD(2)) - YT(L+IADD(1)))*DELTA + YT(L+IADD(1))
	IY(NCROSS) = Y
 25	CONTINUE
	IF(NCROSS .EQ. 0) GO TO 54
	IF(NCROSS .EQ. 2) GO TO 26
	WRITE(IO6, 120) NCROSS, I, J, RLVL, ((RT(II),XT(II),YT(II)),II=
	1  1, 3)
 120	FORMAT('0ERROR>>>>>>>>>>NCROSS = ', I1/' I = ', I3, ', J = ',
	1  I3, ', RLVL = ', F9.4/'0RT', 8X, 'XT', 8X, 'YT'/
	1  (' ', 3(F9.4, 1X), 1X))
	GO TO 54
26  	CONTINUE
	CALL STRCON_DOTLNV(IX(1), IY(1), IX(2), IY(2),DOT,SP)
	GO TO 54
C      CONTOUR LEVEL INTERSECTS ONE VERTEX
 30	CONTINUE
	DO 31 L = 1, 3
 31	IF(IVINT(L)) IV = L
	GO TO (32, 33, 34), IV
 32	IV1 = 2
	IV2 = 3
	GO TO 35
 33	IV1 = 1
	IV2 = 3
	GO TO 35
 34	IV1 = 1
	IV2 = 2
 35	CONTINUE
	IF(RT(IV2) .GT. RT(IV1)) GO TO 36
	IV3 = IV2
	IV2 = IV1
	IV1 = IV3
 36	CONTINUE
	IF((RLVL .LE. RT(IV1)).OR.(RLVL.GE.RT(IV2))) GO TO 54
	DELTA = (RLVL - RT(IV1))/(RT(IV2) - RT(IV1))
	X = (XT(IV2) - XT(IV1))*DELTA + XT(IV1)
	IX(1) = X
	Y = (YT(IV2) - YT(IV1))*DELTA + YT(IV1)
	IY(1) = Y
	IX(2) = XT(IV)
	IY(2) = YT(IV)
	CALL STRCON_DOTLNV(IX(1), IY(1), IX(2), IY(2),DOT,SP)
	GO TO 54
C      CONTOUR LEVEL INTERSECTS TWO VERTICES
 40	CONTINUE
	II = 1
	DO 41 L = 1, 3
	IF(.NOT.IVINT(L)) GO TO 41
	IX(II) = XT(L)
	IY(II) = YT(L)
	II = II + 1
 41	CONTINUE
	CALL STRCON_DOTLNV(IX(1), IY(1), IX(2), IY(2),DOT,SP)
 54	CONTINUE
 53	CONTINUE
 52	CONTINUE
 51	CONTINUE
	RETURN
	END

     	SUBROUTINE STRCON_MARK(X,Y,SIZE,ISYM)
C
C	THIS ROUTINE PUTS ON GRAPHICS MARKERS.
C
     	PARAMETER CIN= 25.4
	COMMON/OFFSET/XOFF,XFAC,YOFF,YFAC

     	XN = (X-XOFF) / XFAC / CIN
     	YN = (Y-YOFF) / YFAC / CIN
     	CALL IFTEXT_MARK( XN, YN, SIZE, 0.0, ISYM )
     	RETURN
     	END

     	SUBROUTINE STRCON_TEXT(X,Y,TEXT,SIZE)
C
C	THIS ROUTINE HANDLES GRAPHTEXT GENERATION.
C
     	PARAMETER CIN= 25.4
     	CHARACTER TEXT*(*)
	COMMON/OFFSET/XOFF,XFAC,YOFF,YFAC

     	SZ = SIZE/CIN
     	OF = 0.5*SZ
     	XN = (X-XOFF) / XFAC / CIN - OF
     	YN = (Y-YOFF) / YFAC / CIN - OF
     	CALL IFTEXT( XN, YN, TEXT, SZ, SZ, 0.0 )
     	RETURN
     	END

