C	WESTERN MICHIGAN UNIVERSITY
C	PUPLOT.FOR (FILE NAME ON LIBRARY DECTAPE)
C	PUPLOT WAS PROGRAMMED BY R.R. BARR III
C	PUPLOT.FOR IS A PLOTTING PROG. CALLED (BY RUNUUO) FROM
C	 PTPLOT.FOR
C	APLB10 PROGS. USED:  FORGEN, IO
C	INTERNAL SUBR. ROUTINES:  AXIS, PAXIS, NICE
C	INTERNAL FUNCTIONS USED:  NINT
C	EXTERNAL SUBR. USED:  (GENERATED BY PTPLOT):  PVPLOT,
C	 PWPLOT, PXPLOT.  ALL 3 SUBR. ARE IN FILE CALLED PVPLOT.
C	FORWMU PROGS. USED:  DEVCHG, RUNUUO, EXISTS, PRINTS, GES,
C	 DEVCHR, GETPPN, JOBNUM
C	SUBR. PVPLOT CONTINAS USER SUPPLIED EQUATIONS IN 
C	 FORTRAN FORM
C	SUBR. PWPLOT TRANSFERS PARAMETERS FROM PTPLOT TO SUBR. PUPLOT
C	 (PTPLOT COMMUNICATES WITH PUPLOT THRU RUNUUO:  THEREFORE
C	 A FORTRAN SUBR. IS REQUIRED.)  THE PARAMETERS THAT
C	 ARE TRANSFERRED ARE IN PTPLOT MAIN PROG. ST. 269+1.
C	PXPLOT TRANSFERS USER'S TITLE INPUT TO PUPLOT FOR
C	 SAME REASON DESCRIBED JUST ABOVE FOR PWPLOT.
C	ABOVE COMMENTS AND RIGHT ADJUSTED COMMENTS PUT IN BY WG
C---------------SEE MAIN PROG. ST. 83, ST. 73-2
	DIMENSION RERUN(7),KAR(10),Y(10),ARRAY(10,0/100),IFMT(16)
C	FOR FORGEN
	DIMENSION MODE(2),JLEN(2)
	DIMENSION PFMT(5),QFMT(5),YA(4)
	DIMENSION MATRIX(0/60,0/100)
C	FOR SUBR PWPLOT & PVPLOT
	COMMON ITYPE,JTYPE,NE,INDVAR,IVAR,ITRAN
C	FOR SUBR IOB
	COMMON/IOB/LEFBK,IRTBK,IALT,MAXPAG,IPAGE,IPAGCT,IDLG,IRSP,II,
     #OUTDV
	COMMON/IOBLK/IDLG1,IRSP1,IDUM(8)
	INTEGER OUTDV
	DOUBLE PRECISION NAMI,NAMO,DVNAM,DEVNAM
C---------------THE NEXT 3 LINES PROBABLY SHOULD BE DELETED.
C	FOR RUN UUO IN 460,460
	DATA RERUN/'RUN P','TPLOT',"555406030140,'0001,','00004',
     1	"565004020100,0/
	DATA KAR/'A','B','C','D','E','F','G','H','I','J'/
	DATA PFMT/'(1X,','15X,','1HI,','51A1',',1HI)'/
	DATA QFMT/'(1X,','F  ,','1H+,','51A1',',1H+)'/
	DATA NDEVI,NDEVO,NDEVT,ITYCH,ISI,ISD,NVSCL,NHSCL,IDLG,IRSP/
     1	4,6,7,0,10,10,5,5,-1,-4/
	IDLG1=IDLG
	IRSP1=IRSP
	TWOPI=2.*3.141592653
	OPEN(UNIT=NDEVT,DEVICE='DSK',MODE='BINARY',FILE='PLTBIN.TMP',
     #ACCESS='SEQINOUT')
C	CALL DEVCHG('DSK',NDEVT)
	CALL IO(1,NDEVO,DVNAM,IDVO,NAMO,IPROJ,IPROG,IBNK)
	IDVO=OUTDV
	IPAGCT=-1
	MODE(1)='F'
	MODE(2)='F'
104	LINCHR=50
	LDNCHR=30
	KNTMAX=0
	PFMT(4)='51A1'
	QFMT(4)='51A1'
	IF(IDVO.EQ.'TTY')GO TO 112
	LDNCHR=LDNCHR*2
	LINCHR=LINCHR*2
	PFMT(4)='101A1'
	QFMT(4)='101A1'
	PFMT(1)=5H('*',
	QFMT(1)=PFMT(1)
112	CALL PWPLOT
	MIS=0
	DO 108 I=0,LDNCHR
	DO 108 J=0,LINCHR
108	MATRIX(I,J)=0
C	TO:   DX  FX  DP  FP
	GO TO(116,85,604,502),(JTYPE-1)*2+ITYPE
C
C---------------DATA PLOTTING WITH RECTANGULAR COORDS.
C
116	CALL IO(0,NDEVI,DEVNAM,IDVI,NAMI,IPROJ,IPROG,IBNK)
	IDVO=OUTDV
	IPAGCT=-1
	CALL FORGEN(IFMT,16,MODE,JLEN,2,2,ISTD,IERR)
	IF(ISTD.EQ.1)IFMT(1)='(2F)'
	WRITE(IDLG,204)
204	FORMAT(' ENTER VERTICAL LIMITS IF DESIRED, ELSE <RETURN>',/)
	READ(IRSP,208)YA(1),YA(3)
208	FORMAT(2F)
	WRITE(IDLG,212)
212	FORMAT(' ENTER HORIZONTAL LIMITS IF DESIRED, ELSE <RETURN>',/)
	READ(IRSP,208)YA(2),YA(4)
	REWIND (NDEVT)
	KI=IVAR
	KD=1
	IF(KI.EQ.1)KD=2
	IF(YA(1).EQ.YA(3))GO TO 214
	XMINN=YA(KI)
	XMAXN=YA(KI+2)
214	IF(YA(2).EQ.YA(4))GO TO 216
	YMINN=YA(KD)
	YMAXN=YA(KD+2)
	REWIND (NDEVT)
216	IF(NDEVI.EQ.NDEVT)GO TO 233
	IF(IDVI.EQ.'TTY')WRITE(IDLG,228)
228	FORMAT(' ENTER DATA',/)
	IF(IDVI.NE.'TTY')WRITE(IDLG,232)
232	FORMAT(' DATA IS BEING READ',/)
C	NDEVI = NDEVT IF THE USER TYPES 'SAME' TO TTY INPUT.
	READ(NDEVI,IFMT,END=670)Y(1),Y(2)
	GO TO 234
233	READ(NDEVT,END=670)Y(1),Y(2)
	GO TO 235
234	IF(IDVI.EQ.'TTY')WRITE(NDEVT)Y(1),Y(2)
235	IF(ITRAN.EQ.'YES')CALL PVPLOT(X,Y)
	XMN=Y(KI)
	XMX=Y(KI)
	YMN=Y(KD)
	YMX=Y(KD)
236	IF(NDEVI.EQ.NDEVT)GO TO 237
	READ(NDEVI,IFMT,END=240)Y(1),Y(2)
	GO TO 238
237	READ(NDEVT,END=240)Y(1),Y(2)
	GO TO 239
238	IF(IDVI.EQ.'TTY')WRITE(NDEVT)Y(1),Y(2)
C---------------FOR DATA PLOTTING BOTH ARGS. ARE INPUT.
C--------------- Y IS MODIFIED AND RETURNED.
239	IF(ITRAN.EQ.'YES')CALL PVPLOT(X,Y)
	IF(Y(KI).LT.XMN)XMN=Y(KI)
	IF(Y(KI).GT.XMX)XMX=Y(KI)
	IF(Y(KD).LT.YMN)YMN=Y(KD)
	IF(Y(KD).GT.YMX)YMX=Y(KD)
	GO TO 236
240	IF(IDVI.EQ.'TTY')REWIND NDEVT
	IF(IDVI.EQ.'TTY')WRITE(IDLG,244)
244	FORMAT(' DATA BEING PROCESSED',/)
	REWIND NDEVI
	IF(YA(KI).EQ.YA(KI+2))CALL NICE(XMN,XMX,ISI,XMINN,XMAXN)
	IF(YA(KD).EQ.YA(KD+2))CALL NICE(YMN,YMX,ISD,YMINN,YMAXN)
	SCALEX=(XMAXN-XMINN)/LINCHR
	SCALEY=(YMAXN-YMINN)/LDNCHR
248	IF(IDVI.EQ.'TTY')GO TO 252
	READ(NDEVI,IFMT,END=260)Y(1),Y(2)
	GO TO 256
252	READ(NDEVT,END=260)Y(1),Y(2)
256	I=(Y(KD)-YMINN)/SCALEY+.5
	J=(Y(KI)-XMINN)/SCALEX+.5
	IF(I.LT.0.OR.J.LT.0.OR.I.GT.LDNCHR.OR.J.GT.LINCHR)GO TO 258
	MATRIX(I,J)=MATRIX(I,J)+1
	IF(MATRIX(I,J).GT.KNTMAX)KNTMAX=MATRIX(I,J)
	GO TO 248
258	MIS=MIS+1
	GO TO 248
260	SCALEK=KNTMAX/FLOAT(MIN0(KNTMAX,10))
	DO 276 I=0,LDNCHR
	DO 276 J=0,LINCHR
	IF(MATRIX(I,J).EQ.0)GO TO 276
	IF(KNTMAX.GT.10)GO TO 264
	K=MATRIX(I,J)
	GO TO 272
264	DO 268 K=1,10
268	IF(SCALEK*K+.5.GE.MATRIX(I,J))GO TO 272
272	MATRIX(I,J)=KAR(K)
276	CONTINUE
	GO TO 86
C
C---------------FUNCTION PLOTTING WITH RECTANGULAR COORDS.
C
85	WRITE(IDLG,304)
304	FORMAT(' ENTER LIMITS FOR DEPENDENT VARIABLE IF DESIRED,',
     1	' ELSE <RETURN>',/)
	READ(IRSP,208)RDA,RDB
308	WRITE(IDLG,312)
312	FORMAT(' ENTER LIMITS FOR INDEPENDENT VARIABLE',/)
	READ(IRSP,208)RIA,RIB
	IF(RIA.NE.RIB)GO TO 318
	WRITE(IDLG,316)
316	FORMAT(' ?RANGE MUST BE NON-ZERO',/)
	GO TO 308
318	WRITE(IDLG,319)
319	FORMAT(' EQUATIONS BEING PROCESSED',//)
C****AM,7.1.1-2,WG,16-DEC-77
	DO 1001 I=1,10
1001	Y(I)=0.
C****END,MAIN PROG.PUPLOT.FOR,STAT.319+4
	SN=SIGN(1.,RIB-RIA)
	SNF=0
	CALL NICE(RIA,RIB,ISI,XMINN,XMAXN)
C---------------FOR FUNCTION PLOTTING X IS INPUT AND Y IS RETURNED.
	CALL PVPLOT(XMINN,Y)
	IF(RDA.EQ.RDB)GO TO 822
	YMINN=RDA
	YMAXN=RDB
822	YMN=Y(1)
	YMX=Y(1)
84	DO 76 J=0,LINCHR
	X=(XMAXN-XMINN)*J/LINCHR+XMINN
	IF(SN.GT.0.AND.(X.LT.RIA.OR.X.GT.RIB))GO TO 76
	IF(SN.LT.0.AND.(X.GT.RIA.OR.X.LT.RIB))GO TO 76
	IF(SNF.EQ.1)GO TO 824
	CALL PVPLOT(X,Y)
	ARRAY(1,0)=Y(1)
	IF(NE.EQ.1)GO TO 823
	DO 83 K=2,NE
	IF(Y(K).LT.YMN)YMN=Y(K)
	IF(Y(K).GT.YMX)YMX=Y(K)
83	ARRAY(K,0)=Y(K)
823	SNF=1
	GO TO 76
824	CALL PVPLOT(X,Y)
	DO 75 K=1,NE
	IF(Y(K).LT.YMN)YMN=Y(K)
	IF(Y(K).GT.YMX)YMX=Y(K)
75	ARRAY(K,J)=Y(K)
76	CONTINUE
	IF(RDA.EQ.RDB)CALL NICE(YMN,YMX,ISD,YMINN,YMAXN)
69	FORMAT(1X,F)
	SINE=RIB-RIA
	XXX=(XMAXN-XMINN)/LINCHR
	DO 74 J=0,LINCHR
	XXY=XXX*J+XMINN
	IF(SINE.GT.0.AND.(XXY.GT.RIB.OR.XXY.LT.RIA))GO TO 74
	IF(SINE.LT.0.AND.(XXY.LT.RIB.OR.XXY.GT.RIA))GO TO 74
	DO 72 K=1,NE
	I=((ARRAY(K,J)-YMINN)/(YMAXN-YMINN))*LDNCHR+.5
C**AM		7.1.1-1, MSL, 13-OCT-77
	IF(I.LT.0.OR.I.GT.LDNCHR)GO TO 73
C**END		PUPLOT MAINLINE, 73 - 2
	MATRIX(I,J)=KAR(K)
	GO TO 72
73	MIS=MIS+1
72	CONTINUE
74	CONTINUE
	XMN=XMINN
	XMX=XMAXN
C
C---------------RECTANGULAR COORDS WITH BOTH DATA AND FUNCTION
C
86	CALL AXIS(LINCHR/50,XMAXN,XMINN,YMAXN,YMINN,MATRIX)
	XS=(XMAXN-XMINN)/NHSCL
	XS1=XS+XMINN
	XS2=2*XS+XMINN
	XS3=3*XS+XMINN
	XS4=4*XS+XMINN
	WRITE(NDEVO,77)
77	FORMAT('1')
	IF(IDVO.EQ.'TTY')WRITE(NDEVO,70)XMINN,XS2,XS4,XS1,XS3,XMAXN
70	FORMAT(8X,F,5X,F,5X,F,/,17X,':',F,4X,':',F,4X,':',F,/,
     1	17X,':',5(9(' '),':'),/,16X,'O+',5(9('-'),'+'),'O')
	IF(IDVO.NE.'TTY')WRITE(NDEVO,170)XMINN,XS2,XS4,XS1,XS3,XMAXN
170	FORMAT('*',7X,F,25X,F,25X,F,/,17X,':',10X,F,14X,':',10X,F,14X,
     1	':',10X,F,/,17X,':',5(19(' '),':'),/,
     1	16X,'O+',5(19('-'),'+'),'O')
	DO 79 I=LDNCHR,0,-1
	IF((I/(LDNCHR/NVSCL)*(LDNCHR/NVSCL)).NE.I)GO TO 180
	SCALE=(YMAXN-YMINN)*I/LDNCHR+YMINN
	WRITE(NDEVO,QFMT)SCALE,(MATRIX(I,J),J=0,LINCHR)
	GO TO 79
180	WRITE(NDEVO,PFMT)(MATRIX(I,J),J=0,LINCHR)
71	FORMAT(' ',15X,'I',51A1,'I')
79	CONTINUE
91	FORMAT(' ',F,'+',51A1,'+')
	IF(IDVO.EQ.'TTY')WRITE(NDEVO,78)XS1,XS3,XMAXN,XMINN,XS2,XS4
78	FORMAT(16X,'O+',5(9('-'),'+'),'O',/,
     1	17X,':',5(9(' '),':'),/,
     1	17X,':',F,4X,':',F,4X,':',F,/,8X,F,5X,F,5X,F)
	IF(IDVO.NE.'TTY')WRITE(NDEVO,178)XS1,XS3,XMAXN,XMINN,XS2,XS4
178	FORMAT('*',15X,'O+',5(19('-'),'+'),'O',/,
     1	17X,':',5(19(' '),':'),/,
     1	17X,':',10X,F,14X,':',10X,F,14X,':',10X,F,/,7X,F,25X,F,25X,F)
	CALL PXPLOT(NDEVO)
	IF(ITYPE.EQ.2)GO TO 991
179	DO 182 I=1,MIN0(KNTMAX,10)
	IF(KNTMAX.GT.10)GO TO 186
	MATRIX(I,0)=I
	MATRIX(I,1)=I
	GO TO 182
186	MATRIX(I,0)=KNTMAX*(I-1)/10.+1.5
	MATRIX(I,1)=KNTMAX*I/10.+.5
182	CONTINUE
	WRITE(NDEVO,184)(KAR(I),MATRIX(I,0),MATRIX(I,1),
     1	I=1,MIN0(KNTMAX,10))
184	FORMAT(/,' SYMBOL    FREQUENCY RANGE OF OCCURRENCE',/,
     1	10(4X,A1,8X,'(',I5,',',I5,')',/))
	IF(KI.EQ.1)KIX='1'
	IF(KI.EQ.2)KIX='2'
	IF(KD.EQ.1)KDX='1'
	IF(KD.EQ.2)KDX='2'
	IF(IVAR.EQ.0)KIX=INDVAR
	IF(IVAR.EQ.0.AND.INDVAR.EQ.'X')KDX='Y'
	IF(IVAR.EQ.0.AND.INDVAR.EQ.'Y')KDX='X'
	IF(JTYPE.EQ.1)WRITE(NDEVO,992)KIX,XMN,XMX,KDX,YMN,YMX
992	FORMAT(36X,'ACTUAL DATA RANGE',/,
     1	' HORIZONTAL VARIABLE(',A1,')',2X,F,1X,F,/,
     1	' VERTICAL VARIABLE(',A1,')',4X,F,1X,F,/)
	WRITE(NDEVO,190)IFMT,IDVI,NAMI
190	FORMAT(' FORMAT:',/,1X,16A5,/,' INPUT - ',A5,':',2A5)
C	EVERYONE ENDS UP HERE
	IF(MIS.NE.0)WRITE(NDEVO,999)MIS
999	FORMAT(/,' NUMBER OF MISSED POINTS = ',I)
991	WRITE(NDEVO,995)
995	FORMAT(///)
	CALL RELEAS(IDLG)
	IF(IDVO.NE.'TTY')WRITE(IDLG,192)
192	FORMAT(/,' PLOT COMPLETE',/)
	IF(ITYPE.EQ.1)GO TO 104
	WRITE(IDLG,68)
68	FORMAT(//,' DO YOU WISH TO ENTER NEW EQUATIONS?(YES OR NO) ',$)
	ACCEPT 67,ANS
67	FORMAT(A3)
	IF(ANS.NE.'YES')GO TO 998
	CALL RELEAS(NDEVO)
	IF(IDVO.EQ.'LPT')CALL PRINTS(NAMO,2,1,1)
C---------------PTPLOT.EXE IS IN SYS:
888	CALL RUNUUO('R PTPLOT ')
C	RUNUUO FOR DISTRIBUTION
C	CALL RUNUUO('RU PTPLOT ')
C	RUNUUO NEVER 'RETURN'S
998	WRITE(IDLG,997)
997	FORMAT(' DO YOU WISH TO PLOT THE SAME FUNCTION(S) WITH',/,
     1	' DIFFERENT LIMITS?(YES OR NO) ',$)
	READ(IRSP,67)ANS
	IF(ANS.EQ.'YES')GO TO 104
	CALL RELEAS(NDEVO)
	IF(IDVO.EQ.'LPT')CALL PRINTS(NAMO,2,1,1)
	CALL EXIT
C
C---------------FUNCTION PLOTTING WITH POLAR COORDS.
C
502	IF(INDVAR.EQ.'T')GO TO 961
964	WRITE(IDLG,962)
962	FORMAT(' ENTER RADIAL LIMIT',/)
	READ(IRSP,972)RLIM
	IF(RLIM.GT.0)GO TO 548
	WRITE(IDLG,963)
963	FORMAT(' ?RADIAL LIMIT MUST BE GREATER THAN ZERO',/)
	GO TO 964
961	WRITE(IDLG,973)
973	FORMAT(' ENTER RADIAL LIMIT IF DESIRED, ELSE <RETURN>',/)
	READ(IRSP,972)RLIM
972	FORMAT(F)
	IF(RLIM.GE.0)GO TO 548
	WRITE(IDLG,963)
	GO TO 961
C	BYPASS FOR NOW
548	GO TO 564
	WRITE(IDLG,974)
974	FORMAT(' ENTER QUADRANT DESIRED, <RETURN> FOR ALL',/)
	READ(IRSP,5561)IQD
5561	FORMAT(I)
	IF(IQD.GE.0.AND.IQD.LE.4)GO TO 564
	WRITE(IDLG,560)
560	FORMAT(' ?RESPONSE ERROR',/)
	GO TO 548
564	XYMX=RLIM
	WRITE(IDLG,319)
C****AM,7.1.1-2,WG,16-DEC-77
	DO 1000 I=1,10
1000	Y(I)=0.
C****END,MAIN PROG. PUPLOT.FOR,STAT.564+5
	IF(INDVAR.EQ.'R')GO TO 520
C
C	THETA INDEPENDENT
C
	IF(RLIM.NE.0)GO TO 966
509	DO 508 J=0,359
	X=TWOPI*J/360.
	CALL PVPLOT(X,Y)
	DO 506 K=1,NE
	IF(ABS(Y(K)*COS(X)).GT.XYMX)XYMX=ABS(Y(K)*COS(X))
	IF(ABS(Y(K)*SIN(X)).GT.XYMX)XYMX=ABS(Y(K)*SIN(X))
506	CONTINUE
508	CONTINUE
966	CALL NICE(0,XYMX,ISD,XYMINN,XYMAXN)
	DO 512 L=0,359
	X=TWOPI*L/360.
	CALL PVPLOT(X,Y)
	DO 510 K=1,NE
	I=(((Y(K)*SIN(X))/XYMAXN)*LDNCHR)/2+.5+LDNCHR/2
	J=(((Y(K)*COS(X))/XYMAXN)*LINCHR)/2+.5+LINCHR/2
	IF(I.GT.LDNCHR.OR.J.GT.LINCHR.OR.I.LT.0.OR.J.LT.0)GO TO 5061
	MATRIX(I,J)=KAR(K)
	GO TO 510
5061	MIS=MIS+1
510	CONTINUE
512	CONTINUE
	GO TO 540
C
C	RADIUS INDEPENDENT
C
520	XYMAXN=XYMX
	DO 536 L=0,LINCHR/2
	X=XYMAXN*L/(LINCHR/2.)
	CALL PVPLOT(X,Y)
	DO 532 K=1,NE
	I=X*SIN(Y(K))/XYMAXN*LDNCHR/2+.5+LDNCHR/2.
	J=X*COS(Y(K))/XYMAXN*LINCHR/2+.5+LINCHR/2
	IF(I.GT.LDNCHR.OR.J.GT.LINCHR.OR.I.LT.0.OR.J.LT.0)GO TO 532
	MATRIX(I,J)=KAR(K)
532	CONTINUE
536	CONTINUE
	GO TO 540
C
C---------------DATA PLOTTING WITH POLAR COORDS.
C
604	CALL IO(0,NDEVI,DEVNAM,IDVI,NAMI,IPROJ,IPROG,IBNK)
	IDVO=OUTDV
	IPAGCT=-1
	CALL FORGEN(IFMT,16,MODE,JLEN,2,2,ISTD,IERR)
	IF(ISTD.EQ.1)IFMT(1)='(2F)'
	REWIND (NDEVT)
	WRITE(IDLG,973)
C	RLIM IF DESIRED
	READ(IRSP,972)RLIM
C	BYPASS FOR NOW
608	GO TO 612
	WRITE(IDLG,974)
	READ(IRSP,5561)IQD
	IF(IQD.GE.0.AND.IQD.LE.4)GO TO 612
	WRITE(IDLG,560)
	GO TO 608
612	KR=IVAR
	KT=2
	IF(KR.EQ.2)KT=1
	IF(IDVI.EQ.'TTY')WRITE(IDLG,228)
	IF(IDVI.NE.'TTY')WRITE(IDLG,232)
	REWIND(NDEVT)
C	NDEVI=NDEVT IF USER TYPES 'SAME' TO TTY INPUT
	IF(NDEVI.EQ.NDEVT)GO TO 614
	READ(NDEVI,IFMT,END=670)Y(1),Y(2)
	GO TO 616
614	READ(NDEVT,END=670)Y(1),Y(2)
	GO TO 618
616	IF(IDVI.EQ.'TTY')WRITE(NDEVT)Y(1),Y(2)
618	IF(RLIM.NE.0)GO TO 620
	IF(ITRAN.EQ.'YES')CALL PVPLOT(X,Y)
	TTMP=AMOD(Y(KT),TWOPI)
	IF(TTMP.LT.0)TTMP=TWOPI+TTMP
	TMN=TTMP
	TMX=TTMP
	RMN=Y(KR)
	RMX=Y(KR)
620	IF(NDEVI.EQ.NDEVT)GO TO 621
	READ(NDEVI,IFMT,END=624)Y(1),Y(2)
	GO TO 622
621	READ(NDEVT,END=624)Y(1),Y(2)
	GO TO 623
622	IF(IDVI.EQ.'TTY')WRITE(NDEVT)Y(1),Y(2)
623	IF(RLIM.NE.0)GO TO 620
	IF(ITRAN.EQ.'YES')CALL PVPLOT(X,Y)
	TTMP=AMOD(Y(KT),TWOPI)
	IF(TTMP.LT.0)TTMP=TWOPI+TTMP
	IF(TTMP.LT.TMN)TMN=TTMP
	IF(TTMP.GT.TMX)TMX=TTMP
	IF(Y(KR).LT.RMN)RMN=Y(KR)
	IF(Y(KR).GT.RMX)RMX=Y(KR)
	GO TO 620
624	IF(IDVI.EQ.'TTY')REWIND NDEVT
	IF(IDVI.EQ.'TTY')WRITE(IDLG,244)
	REWIND NDEVI
	RMAXN=RLIM
	IF(RLIM.EQ.0)CALL NICE(RMN,RMX,ISD,RMINN,RMAXN)
628	IF(IDVI.EQ.'TTY')GO TO 632
	READ(NDEVI,IFMT,END=640)Y(1),Y(2)
	GO TO 636
632	READ(NDEVT,END=640),Y(1),Y(2)
636	IF(ITRAN.EQ.'YES')CALL PVPLOT(X,Y)
	I=Y(KR)*SIN(Y(KT))/RMAXN*LDNCHR/2+.5+LDNCHR/2
	J=Y(KR)*COS(Y(KT))/RMAXN*LINCHR/2+.5+LINCHR/2
	IF(I.GT.LDNCHR.OR.J.GT.LINCHR.OR.I.LT.0.OR.J.LT.0)GO TO 638
	MATRIX(I,J)=MATRIX(I,J)+1
	IF(MATRIX(I,J).GT.KNTMAX)KNTMAX=MATRIX(I,J)
	GO TO 628
638	MIS=MIS+1
	GO TO 628
640	SCALEK=KNTMAX/FLOAT(MIN0(KNTMAX,10))
	DO 660 I=0,LDNCHR
	DO 660 J=0,LINCHR
	IF(MATRIX(I,J).GT.10)GO TO 644
	K=MATRIX(I,J)
	GO TO 656
644	DO 648 K=1,10
648	IF(SCALEK*K+.5.GE.MATRIX(I,J))GO TO 656
656	MATRIX(I,J)=KAR(K)
660	CONTINUE
C
C	ALL POLAR PLOTTING METHODS PRINT HERE
C
540	CALL PAXIS(LINCHR/50,MATRIX)
	FCH=' '
	IF(IDVO.NE.'TTY')FCH='*'
	BLANK=' '
	MARGIN=15
	WRITE(NDEVO,77)
	IF(IDVO.EQ.'TTY')WRITE(NDEVO,558)XYMAXN
	IF(IDVO.NE.'TTY')WRITE(NDEVO,5591)XYMAXN
558	FORMAT(33X,F)
5591	FORMAT(58X,F)
	JMAX=LINCHR
	DO 971 I=LDNCHR,0,-1
	IF(IDVO.NE.'TTY')GO TO 554
	DO 552 J=LINCHR,1,-1
	JMAX=J
	IF(MATRIX(I,J).NE.' ')GO TO 554
552	CONTINUE
554	IF(I.EQ.LDNCHR/2)GO TO 556
	WRITE(NDEVO,544)FCH,(BLANK,J=1,MARGIN),(MATRIX(I,J),J=0,JMAX)
544	FORMAT(132A1)
	GO TO 971
556	WRITE(NDEVO,559)FCH,XYMAXN,(MATRIX(I,J),J=0,JMAX)
559	FORMAT(A1,F,101A1)
971	CONTINUE
	IF(IDVO.EQ.'TTY')WRITE(NDEVO,558)XYMAXN
	IF(IDVO.NE.'TTY')WRITE(NDEVO,5591)XYMAXN
	CALL PXPLOT(NDEVO)
	IF(IQD.NE.0)GO TO 952
	IQS=0
	IQD=359
952	IF(RMAXN.NE.0)XYMAXN=RMAXN
	TICD=XYMAXN/2.5/(LINCHR/50)
	IF(ITYPE.EQ.2)WRITE(NDEVO,954)INDVAR
954	FORMAT(' THE INDEPENDENT VARIABLE IS ',A1,/)
	WRITE(NDEVO,951)TICD
951	FORMAT(' DISTANCE BETWEEN "+" SIGNS IS ',F)
	IF(IQS.NE.0.OR.IQD.NE.359)WRITE(NDEVO,953)IQS,IQD
953	FORMAT(' RANGE OF THETA IS ',I3,',',I3)
	IQD=0
	IF(ITYPE.EQ.1)GO TO 179
	GO TO 991
670	WRITE(IDLG,672)
672	FORMAT('?NO DATA IN THIS FILE',/)
	GO TO 104
	END
C---------------PUTS IN POLAR AXES MARKS WHICH DO NOT CONFLICT
C--------------- WITH DATA ALREADY PRESENT IN MATRIX
C---------------N, MATRIX ARE INPUT.  MATRIX IS MODIFIED.
	SUBROUTINE PAXIS(N,MATRIX)
	DIMENSION MATRIX(0/60,0/100),IPTS(7,2)
	DATA ((IPTS(I,J),J=1,2),I=1,7)
     1	/4,7,9,14,13,22,17,28,21,35,25,42,30,49/
	DO 1 I=0,5*N/2
	DO 2 J=-1,1,2
	IF(MATRIX(N*15*J*6,25*N).EQ.0)MATRIX(N*15+I*6*J,25*N)='+'
	IF(MATRIX(N*15,25*N+I*10*J).EQ.0)MATRIX(N*15,25*N+I*10*J)='+'
2	CONTINUE
1	CONTINUE
	DO 7 I=0,30*N
7	IF(MATRIX(I,25*N).EQ.0)MATRIX(I,25*N)='*'
	DO 8 J=0,50*N
8	IF(MATRIX(15*N,J).EQ.0)MATRIX(15*N,J)='*'
5	DO 98 I=1,7
	IF(N.EQ.1.AND.I.GT.3)GO TO 95
	DO 98 IB=-1,1,2
	DO 98 JB=-1,1,2
	IA=IPTS(I,1)*IB
	JA=IPTS(I,2)*JB
	IF(MATRIX(IA+15*N,JA+25*N).EQ.0)MATRIX(IA+15*N,JA+25*N)='+'
98	CONTINUE
95	DO 99 I=0,30*N
	DO 99 J=0,50*N
	IF(MATRIX(I,J).EQ.0)MATRIX(I,J)=' '
99	CONTINUE
	RETURN
	END
C---------------PUTS IN X AND Y AXES MARKS WHICH DO NOT CONFLICT
C--------------- WITH DATA ALREADY PRESENT IN MATRIX.  ALL ARGS.
C--------------- ARE MODIFIED.  MATRIX IS MODIFIED.
	SUBROUTINE AXIS(MULT,XMAXN,XMINN,YMAXN,YMINN,MATRIX)
	DIMENSION MATRIX(0/60,0/100)
	IF(XMAXN.LT.0.OR.XMINN.GT.0)GO TO 402
	K=ABS(XMINN/(XMAXN-XMINN))*50*MULT+.5
	DO 401 I=0,30*MULT
	IF(MATRIX(I,K).EQ.0)MATRIX(I,K)='*'
401	CONTINUE
402	IF(YMAXN.LT.0.OR.YMINN.GT.0)GO TO 404
	K=ABS(YMINN/(YMAXN-YMINN))*30*MULT+.5
	DO 403 J=0,50*MULT
	IF(MATRIX(K,J).EQ.0)MATRIX(K,J)='*'
403	CONTINUE
404	DO 405 I=0,30*MULT
	DO 405 J=0,50*MULT
	IF(MATRIX(I,J).EQ.0)MATRIX(I,J)=' '
405	CONTINUE
	RETURN
	END
C	THIS FUNCTION IS USED BY NICE TO PROVIDE AN INTEGER
C	FUNCTION THAT ALLWAYS ROUNDS NEGITIVE
C	I.E.
C		2=NINT(2.5)
C		-3=NINT(-2.5)
C
C---------------X IS INPUT
	FUNCTION NINT(X)
	NINT=X
	IF(FLOAT(NINT).GT.X)NINT=NINT-1
	RETURN
	END
C
C	NICE - A ROUTINE TO SET (FROM A HUMAN POINT OF VIEW)
C	INTERVALS FOR A PLOT SCALE.
C
C	WRITTEN BY RUSSELL R. BARR III - WMU COMPUTER CENTER- KALAMAZOO
C	DATE: MARCH 1973
C
C	THIS COPY OF NICE IS MADE TO USE THE METHOD EXPLAINED IN
C	'DIGITAL TIME SERIES ANALYSIS', BY OTNES AND ENOCHSON
C	QA 280.087(PYS SCI) PG. 65
C	IT BEARS IMPROVEMENTS AND CORRECTIONS.
C---------------XMIN, XMAX, S ARE INPUT.  E, F ARE OUTPUT.
	SUBROUTINE NICE(XMIN,XMAX,S,E,F)
	INTEGER S
	DIMENSION DI(10)
	DATA N,DI/5,1.5,2.,2.5,5.,10.,5*0./
	IFLG=0
	IF(XMAX.NE.XMIN)GO TO 6
	E=XMIN
	F=XMAX
	RETURN
6	IF(XMAX.GE.XMIN)GO TO 4
	X=XMAX
	XMAX=XMIN
	XMIN=X
	IFLG=1
4	XX=(XMAX-XMIN)/S
	X=ALOG10(XX)
	I=NINT(X)
	C=XX/10.**I
D	TYPE 99,XX,I,C
99	FORMAT(' XX1='F,' I='I,' C='F)
	C=(C+100.)-100.
	DO 100 J=N,1,-1
	IF(C.LE.DI(J))D=DI(J)
D	TYPE 93,C,D,DI(J),J
93	FORMAT(3F,I2)
D	TYPE 94,C,DI(J)
94	FORMAT(2O)
100	CONTINUE
D	TYPE 98,D
1	XX=D*10.**I
	J=NINT(XMIN/XX)
	E=XX*J
	F=E+XX*S
D	TYPE 98,D,XX,J,E,F
98	FORMAT(' D='F,' XX2='F,' J='I,' E='F,' F='F)
	IF(XMAX.LE.(F+100.)-100.)GO TO 5
	DO 110 J=1,N-1
	IF(D.NE.DI(J))GO TO 110
	D=DI(J+1)
D	TYPE 98,D
	GO TO 1
110	CONTINUE
	D=DI(1)
	I=I+1
D	TYPE 97,D,I
97	FORMAT(' D='F,' I='I)
	GO TO 1
5	IF(IFLG.EQ.0)RETURN
	X=XMAX
	XMAX=XMIN
	XMIN=X
	X=F
	F=E
	E=X
	RETURN
	END

