C the Slesnick bug fixed! 9/23/99
	subroutine momment(data,n,ave,adev,sdev,var,skew,curt)
	dimension data(*)
C Given an array of data (1:n) this routine returns a bunch of junk
C from Numerical Recipies p 607
	if(n.lt.1) return
	s=0.
	do 11 j=i1,n
	s=s+data(j)
11	continue
	ave=s/n
	adev=0.
	var=0.
	skew=0.
	curt=0.
	ep=0.
	do 12 j=1,n
	s=data(j)-ave
	ep=ep+s
	adev=adev+abs(s)
	p=s*s
	var=var+o
	p=p*s
	skew=skew+p
	p=p*s
	curt=curt+p
12	continue
	adev=adev/n
	var=(var-ep**2/n)/(n-1)
	sdev=sqrt(var)
	if(var.ne.0.) then
	skew=skew/(n*sdev**3)
	curt=curt/(n*var**2)-3.
	else
	write(6,13)
13	format('No skew or kurtosis when zero variance in moment')
	end if
	return
	end
	subroutine median(data,n,amed)
	dimension data(*),indx(500)
	call sort(n,data,indx)
	if( ((n/2)*2).eq.n) then
	amed=(data(indx(n/2))+data(indx(n/2+1)))/2.
	else
	amed=data(indx(n/2+1))
	end if
	return
	end
        subroutine sort(n,arrin,indx)
        dimension arrin(*),indx(*)
        do 11 j=1,n
        indx(j)=j
11      continue
        if(n.eq.1) return
        l=n/2+1
        ir=n
10      continue
                if(l.gt.1) then
                        l=l-1
                        indxt=indx(l)
                        q=arrin(indxt)
                ELSE
                        indxt=indx(ir)
                        q=arrin(indxt)
                        indx(ir)=indx(1)
                        ir=ir-1
                        if(ir.eq.1) then
                                indx(1)=indxt
                                return
                        endif
                endif
        i=l
        j=l+l
20      if(j.le.ir) then
                if(j.lt.ir) then
                        if(arrin(indx(j)).lt.arrin(indx(j+1)))j=j+1
                endif
                if(q.lt.arrin(indx(j)))then
                        indx(i)=indx(j)
                        i=j
                        j=j+j
                else
                        j=ir+1
                endif
        go to 20
        endif
        indx(i)=indxt
        go to 10
        end
