SUBROUTINE HISTOGRM(SAMP,NSAMP,NCL,XLL,XUU,barwidth,IO)
	common /hist/ hd
	common /grades/ jj, gr
	INCLUDE 'SYS$CHELIB:PLTGEN.FOR'
	character*3 legend
	character*1 lg(11)*2,LLG*22
	character*9 datebuf
	character*7 hd(46)
	integer reject
	dimension nn(1000),fn(1000),XM(1001)
	dimension samp(1),xpt(4),ypt(4)
	dimension gr(11)
	equivalence (lg(1),llg)
	data llg/'E D D+C-C C+B-B B+A-A '/
c
	call date(datebuf)
	smin = samp(1)
	smax = samp(1)
	xm(1) = xll
	xm(12) = xuu
c
	do i=1,100
	nn(i)=0
	enddo
c
	reject=0
c
	type *,' Do you want to'
	type *,' 1. Use the ranges given in cutoff values'
	type *,' 2. Choose equal ranges between maximum'
	type *,'    and minimum values entered above'
	type *,' 3. Consider other options'
50	type *,' Enter 1, 2 or 3'
	read (5,*,err=50)icode
	if(icode.lt.1 .or. icode.gt.3)goto 50
	goto (60,70,80),icode
c
60	do j=2,ncl
	xm(j)=gr(j-1)
	enddo
	goto 150
c
70	del=(xuu-xll)/float(ncl)
	do n=1,ncl+1
	xm(n)=xll+del*float(n-1)
	enddo
	goto 150
c
150	do  n=1,nsamp
	if(samp(n).gt.xuu.or.samp(n).lt.xll)then
		reject=reject+1
 		go to 1100
	endif
c
	do kk=1,ncl+1
	if(samp(n).ge.xm(kk) .and. samp(n).lt.(xm(kk+1)+.0001))goto 90
	enddo
90	i=kk
c
	nn(i)=nn(i)+1
		if(samp(n).gt.smax)then
			smax=samp(n)
		  else if(samp(n).lt.smin)then
			smin=samp(n)
		  else
		endif
1100	enddo
c
	fmax=0.
	nnmax=0
	do n=1,ncl+1
	if(nn(i).gt.nnmax)nnmax=nn(i)
	fn(n)=float(nn(n))/float(nsamp)*100.
	if(fn(n).gt.fmax)then
		fmax=fn(n)
		nmax=n
		         else
		continue
	endif
	enddo
c
	n1=1
	inc=6
	do nnn=1,ncl,inc
	n2=n1+inc-1
	if(n2.gt.ncl)n2=ncl
	write(io,11)(xm(I),I=n1,n2+1)
	write(io,12)(NN(I),I=n1,n2)
	write(io,14)(lg(i),i=n1,n2)
        write(io,13)(fn(I),I=n1,n2)      
	n1=n1+inc
	enddo
11	format(// f6.2,11f9.2)
12	format(11I9)
13	format(2x,11f9.2)
14	format(11A9)
c	
	if(reject.ne.0)write(io,20)reject
20	format(////2x,i3,'  samples out of bound for histogram')
c	
	call stat(samp,nsamp,xmean,sdev,io)
c
	fac=0.49
	xorgin=1.5
	yorgin=1.5
	xs1=xorgin
	ys1=yorgin
	call plt$gen(1,2,fac)
	ky=ifix(fmax/10.+1.5)
	xnmax=float(nnmax)
	xfac=(xuu-xll)/xaxlen
	yfac=10.*ky/yaxlen
c
	CALL plt$ax(xll, xuu, ncl/2 ,1,  
     *                                  0., 10.*ky, ky, 0, 0,0,1)
	call plt$tl('Frequency Distribution',hd(jj),'Percent',1)
	ypt(1)=0.0
	ypt(4)=0.0
c
	do n=1,ncl
	xpt(1)=xm(n)
	xpt(2)=xm(n)
	xpt(3)=xm(n+1)
	xpt(4)=xm(n+1)
	ypt(2)=fn(n)
	ypt(3)=fn(n)
	call plt$pt(xpt,ypt,4,0,0,0,1)
c
	xnum = float(nn(n))
	inum=xnum
c
	if(inum.gt.0)then
	legend=lg(n)//'('
	xlg=(xpt(1)+xpt(3))/2.0
	ylg=(ypt(3)+0.1*yfac)
	call plt$lg(xlg,ylg,0,legend,90.,1)
	YLG=YLG+YFAC*3.*.36
	call plt$num(xlg,ylg,xnum,0,90.0,1)
	YLG=YLG+YFAC*2.*.36
	call plt$lg(xlg,ylg,0,')',90.,1)
	endif
c
	nj=nn(n)-1
	do jjj=1,nj
	yloc=fn(n)/nn(n)*float(jjj)
	x1=(xm(n)-xll)/xfac
	x2=(xm(n+1)-xll)/xfac
	call plot(x1,yloc/yfac,3)
	call plot(x2,yloc/yfac,2)
	enddo
c
9000	enddo
	y1=0.1
	X1=(XMEAN-XLL)/XFAC-0.1
	X2=(XMEAN-XLL)/XFAC+0.1
	CALL PLOT(X1,-Y1,3)
	CALL PLOT(X2,Y1,2)
	CALL PLOT(X2,-Y1,3)
	CALL PLOT(X1,Y1,2)
c
	CALL PLT$HD(DATEBUF,1)
	call tsend
	call trans(23,24)
	type *,' Do you want to screendump this file (Y or N)? '
	read(5,1015)an
1015	format(a1)
	if(an.eq.'y'.or.an.eq.'Y')call scrdump
	call vecmod
	call newpag
	call plt$dn
	call trans(1,24)
c
80	RETURN
25	format(1x,11(f6.1,1x))
	end
	subroutine stat(s,n,ave,sd,io)
	dimension s(1)
c	character*16 fl
c	common/name/fl
	SUM=0.
	DO I=1,n
	SUM=SUM+s(i)
	ENDDO
C
	ave=sum/float(n)
	ss=0.
	do i=1,n
	sq=(s(i)-ave)**2
	ss=ss+sq
	enddo
	sd=sqrt(ss/(float(n-1)))
c	write(55,*)n,ave,sd,fl
	write(io,15)n,ave,sd
15	format(//3x,'n = ',i4,5x,'mean = ',f5.2,5x,'s.d. = ',f5.2)
c
	return
	END