	PROGRAM PDPLOT
C
C PLOT OBSERVED AND CALCULATED POWDER DIFFRACTION PROFILE FOR PROF
C
C before ruinning, run a command procedure to get files right.
c $ASSIGN RBC.P2PLO FOR095
c $ASSIGN RBC.HKL   FOR060
c $ASSIGN RBC.P1INP FOR040
c $ASSIGN RBC.PROFIN FOR043
c $ASSIGN RBC.COORD FOR042

C Modified 7/31/91 PWS, to take data from .PROFIN, so that excluded 
C regions are plotted.  Also, take unit cell parameters from refined .COORD

	PARAMETER (MAXDATA=8000)
	DIMENSION IOBS(10),ICAL(10)
        DIMENSION TTHETA(1000), YBUF(MAXDATA), POSBUF(MAXDATA)
	DIMENSION POS(MAXDATA),YCAL(MAXDATA),YOBS(MAXDATA),YDIF(MAXDATA)
	CHARACTER*80 PTITLE,STTIME*8,STDATE*9
	LOGICAL HARDF
	INTEGER LO,LI,ICRT
	LOGICAL ANSI
	DIMENSION XBAR(2),YBAR(2)
C
	NOPLOT = 0
	LO = 1
	LI = 1
	OPEN(UNIT=LO,FILE='SYS$OUTPUT',STATUS='NEW')
	OPEN(UNIT=LI,FILE='SYS$INPUT',STATUS='NEW')
	ICRT = 0
	CALL SETCRT(ICRT,LO,LI,ANSI)
	CALL ERRSET(29,.TRUE.,,.TRUE.)
C
C GET TITLE AND SCAN ANGLES
C
	READ(95,'(A80)') PTITLE
	READ(95,'(F8.3,F8.5,F8.3)') ANGSTR,ANGSTP,ANGEND
C
C GET CELL DIMENSIONS AND WAVELENGTH
C
	READ(40,'(//)')		!throw away 2 lines
C	READ(40,'(6F12.4)')A
      	READ(40,'(F8.4)') FLAME	! x-ray lambda
	CLOSE (40)
c	read lattice constants from COORD file
110	read(42,'()', end=120)
	go to 110
120	do i=1,4
		backspace 42
	enddo
	read (42,*) AL, BL, CL, DL, EL, FL	! REFINED LATTICE PARAMS
	type 1000, AL, BL, CL, DL, EL, FL
1000	format(' LATTICE PARAMS: ',6f10.7)
C
C GET REFLECTION LIST
C
	IM=0
	DO I=1,1000
	  READ(60,'(8X,3I8)',END=10) IH, IK, IL
	  DM2 = IH*IH*AL + IK*IK*BL + IL*IL*CL + IK*IL*DL + IH*IL*EL + IH*IK*FL
	  PAM=0.5*FLAME*SQRT(DM2)
	  IF(PAM.GT.1.0) GOTO 5
          IM=IM+1
  	  TTHETA(IM)=2.0*ASIND(PAM)
5	  CONTINUE
	ENDDO
10	NREF=I-1
C
	N=0
	POSMIN=4.5
 	POSMAX=35.
	YMIN=10000.
	YMAX=0.
	YDMIN=10000.
	YDMAX=0.
	DO K=1,800
	  READ(95,'(10I8)',END=20) IOBS
  	  READ(95,'(10I8)',END=20) ICAL
	  DO I=1,10
	    N=N+1
	    POS(N) =ANGSTR+(N-1)*ANGSTP
	    IF(POS(N).GT.ANGEND) GOTO 20
	    YCAL(N)=ICAL(I)
	    YOBS(N)=IOBS(I)
	    YDIF(N)=IOBS(I)-ICAL(I)
	    POSMIN=MIN(POS(N),POSMIN)
	    POSMAX=MAX(POS(N),POSMAX)
	    YMIN=MIN(YCAL(N),YMIN)
	    YMIN=MIN(YOBS(N),YMIN)
	    YMAX=MAX(YCAL(N),YMAX)
   	    YMAX=MAX(YOBS(N),YMAX)
	    YDMIN=MIN(YDIF(N),YDMIN)
	    YDMAX=MAX(YDIF(N),YDMAX)
       	  ENDDO
	ENDDO
C ABOVE HAS ESTABLISHED LIMITS TO ENCOMPASS DATA AND FIT IN INCLUDED REGIONS

C NOW READ RAW DATA FROM *.PROFIN
20	READ(43, '(/)')		!SKIP FIRST 2 LINES
	DO K=1,800
C	TYPE 1003, K, IFIX(YOBS(10*K - 9))
	  READ(43, '(10F8.0)', END = 21) IOBS
	  DO I=1,10
	    NOBS = 10*K + I - 10
	    IF(NOBS .GT. N) GO TO 21
	    YOBS(NOBS) = IOBS(I)
	  ENDDO
C	TYPE 1003, K, IFIX(YOBS(10*K - 9))
1003	FORMAT(2I10)
	ENDDO
	    
C
21	HARDF = .FALSE.
30	IF (HARDF) THEN
	  CALL OPNHD(LO,LI)
	ELSE
	  CALL BGRAF(ICRT,LO,LI)
	ENDIF
C
C LINES AT POSITIONS OF ALLOWED PEAKS
C
	IF(NOPLOT) GO TO 401
	CALL PGVPORT(0.05,0.95,0.25,0.30)
	CALL PGWINDOW(POSMIN,POSMAX,0.,1.)
	YBAR(1)=0.2
	YBAR(2)=0.8
	DO K=1,IM
	  XBAR(1)=TTHETA(K)
	  XBAR(2)=TTHETA(K)
	  CALL PGLINE(2,XBAR,YBAR)
 	ENDDO
401	CONTINUE
C
C  
	IF(NOPLOT) GO TO 402
	CALL PGVPORT(0.05,0.95,0.3,1.0)
	CALL PGWINDOW(POSMIN,POSMAX,YMIN,YMAX)
	CALL PGBOX('BCST',0.0,0,'BCNST',0.0,0)
C	CALL DATE(STDATE)
C	CALL PGMTEXT('T',-1.5,0.97,1.0,STDATE)
C	CALL PGMTEXT('T',-2.5,0.97,1.0,PTITLE)
	CALL PGPOINT(N,POS,YOBS,1)
402	CONTINUE
C
C  PLOT YCAL ONLY FOR NON-EXCLUDED REGIONS
C  LOOP ON NON-EXCLUDED REGIONS.
	NSTART = 1
	DO WHILE(.TRUE.)
	  DO WHILE(.TRUE.)
	    IF(NSTART .GE. N) GO TO 301		! DONE
	    IF(YCAL(NSTART) .NE. 0) GO TO 302	!START HERE
	    NSTART = NSTART + 1
	  ENDDO
	  GO TO 301
302	  NSTOP = NSTART+1
	  DO WHILE (NSTOP .LE. N .AND. YCAL(NSTOP) .NE. 0)
	    NSTOP = NSTOP + 1
	  ENDDO
	  NSTOP = NSTOP-1
	  NN = 1
	  DO NNN = NSTART, NSTOP
	    YBUF(NN) = YCAL(NNN)
	    POSBUF(NN) = POS(NNN)
	    NN = NN + 1
	  ENDDO
	IF(NOPLOT .EQ. 0) GO TO 404
	TYPE 1001, NSTART, NSTOP
1001	FORMAT(' ID PLOT FROM ', I5, ' TO ', I5)
404	CONTINUE
	  IF(NOPLOT) GO TO 403
	  CALL PGLINE(NSTOP-NSTART+1, POSBUF, YBUF)
403	  CONTINUE
	NSTART = NSTOP+1
	ENDDO
301	CALL PGSLW(1)

C
C  PLOT DIFFERENCES ONLY FOR NON-EXCLUDED REGIONS

	CALL PGVPORT(0.05,0.95,0.10,0.25)
	CALL PGWINDOW(POSMIN,POSMAX,YDMIN,YDMAX)
	CALL PGBOX('BCNST',0.0,0,'BCST',0.0,0)
	NSTART = 1
	DO WHILE(.TRUE.)
	  DO WHILE(.TRUE.)
	    IF(NSTART .GE. N) GO TO 501		! DONE
	    IF(YCAL(NSTART) .NE. 0) GO TO 502	!START HERE
	    NSTART = NSTART + 1
	  ENDDO
	  GO TO 501
502	  NSTOP = NSTART+1
	  DO WHILE (NSTOP .LE. N .AND. YCAL(NSTOP) .NE. 0)
	    NSTOP = NSTOP + 1
	  ENDDO
	  NSTOP = NSTOP-1
	  NN = 1
	  DO NNN = NSTART, NSTOP
	    YBUF(NN) = YOBS(NNN) - YCAL(NNN)
	    POSBUF(NN) = POS(NNN)
	    NN = NN + 1
	  ENDDO
	  IF(NOPLOT) GO TO 405
	  CALL PGPOINT(NSTOP-NSTART+1, POSBUF, YBUF, 1)
405	  CONTINUE
	NSTART = NSTOP+1
	ENDDO
501	CONTINUE

C REPEAT EVERYTHING TO MAKE HARD COPY
	IF (HARDF) THEN
	  HARDF = .FALSE.
	  CALL CLHD(LO,LI)
	ELSE
	  WRITE(*,'(A)') '$DO YOU WANT A HARDCOPY? [N]: '
	  READ (*,'(A)',END=999) IANS
	  CALL EGRAF(ICRT,LO,LI)
	  IF(IANS.EQ.'Y' .or. IANS.EQ.'y') THEN 
	    HARDF = .TRUE.
	    GOTO 30
	  ENDIF
	ENDIF
999	CALL EXIT
	END
