      subroutine qdf_gnu(iug,n,x) ! bernd berg, june 17 2002.
c gnuplot for the peaked cumulative distribution function (qcdf).
c the array x has to be sorted on input. 
      include '../../libs/fortran/implicit.sta'
      include '../../libs/fortran/constants.par'
      character ci*2
      dimension x(n)
      data icnt /0/
      save icnt
c
      icnt=icnt+1
      if(icnt.gt.99) stop 'df_gnu: inct=99 exhausted!'
      write(ci,'(i2.2)') icnt
c
c gnuplot script:
      open(iug,file='qdf'//ci//'.plt',status='unknown',form='formatted')
      write(iug,*) 'set noyzeroaxis'
      write(iug,*) 'set yrange [0:0.525]'
      write(iug,*) 'plot "qdf'//ci//'.d" using 1:3 with l lt 1'
      write(iug,*) 'pause -1'
      close(iug)
c
      if(icnt.gt.1) then ! create qdf.plt to plot all histograms
      open(iug,file='qdf.plt',status='unknown',form='formatted')
      write(iug,*) 'set noyzeroaxis'
      write(iug,*) 'set yrange [0:0.525]'
      write(iug,*) 'plot "qdf01.d" using 1:3 with l lt 1,\\'
      do i=2,(icnt-1)
        write(iug,101) i,i
      end do
101   format('      "qdf',i2.2,'.d" using 1:3 with l lt ',i2,',\\')
      write(iug,103) icnt,icnt
103   format('      "qdf',i2.2,'.d" using 1:3 with l lt ',i2)
      write(iug,*) 'pause -1'
      close(iug)
      endif
c
      open(iug,file='qdf'//ci//'.d',status='unknown',form='formatted')
      write(iug,'(3g16.7)') x(1),zero,zero
      do i=1,(n-1)
        if(x(i+1).lt.x(i)) stop "qdf_gnu: data not sorted."
        f=(i*one)/(n*one)
        fq=min(f,one-f)
        write(iug,'(3g16.7)') x(i),  f,fq
        write(iug,'(3g16.7)') x(i+1),f,fq
      end do
      write(iug,'(3g16.7)') x(n),one,zero
      close(iug)
c
      return
      end
