envelop.f

来自「这是一个希尔伯特变换的程序包」· F 代码 · 共 34 行

      subroutine hilbert(lx,x)
      integer   jsz, item, ierr, iabort
      dimension h(129),x(*)
      real  y
      pointer (wky, y(1))
      data  ierr/0/
      data  iabort/0/

      lh=129
      item = lx+lh+1
      call sizefloat(jsz)
      call galloc (wky, jsz*item, ierr, iabort)
      if (ierr .ne. 0) then
         call ppc('cross (hilbert): FATAL ERROR')
         call ppc('Unable to allocate bytes',jsz*item)
         return
      endif
      lh2=(lh-1)/2
      call vclr(h,1,lh)
      do  10  i = 1, lh2, 2
         d=i
         h(lh2+1-i) = -1./d
         h(lh2+1+i) =  1./d
   10 continue
      call dotpr (h,1,h,1,s,lh)
      s=1./sqrt(s)
      do  11  i = 1,lh
         h(i) = s * h(i)
   11 continue
      call fold(lh,h,lx,x,ly,y)
      call vmov(y(lh2+1),1,x,1,lx)
      return
      end