      subroutine knotins(tau,t,k,j,left,alpha)
c
c     This subroutine is an implementation of algorithm 1
c     for calculating all non-vanishing discrete B-splines
c     according to a recurrence relation by Cohen,Lyche,Riesenfeld:
c     Discrete B-splines and subdivision Techniques in
c     computer-aided geometric design and computer graphics,
c     Computer graphics and image processing 14, 87-111 (1980)
c
c     INPUT:
c           tau    = initial sequence of knots
c           t      = knot refinement, may include some of the tau knots
c           k      = order of the spline space
c           j,left = integers such that tau(left).le.t(j).lt.tau(left+1)
c     OUTPUT:
c           alpha  = vector of k non-vanishing discrete B-splines
c     WARNING:
c           division by zero ocurrs if the condition
c           tau(left).le.t(j).lt.tau(left+1) is not satisfied
c
      implicit double precision (a-h,o-z)
      dimension alpha(left-k+1:left),tau(1),t(1)
      alpha(left)=1.d0
      l2=left
      do 1 ir=1,k-1
         beta_1=0.d0
         tj=t(j+ir)
         do 2 i=l2,left
            d1=tj-tau(i)
            d2=tau(i+ir)-tj
            beta=alpha(i)/(d1+d2)
            alpha(i-1)=d2*beta+beta_1
            beta_1=d1*beta
2        continue
         alpha(left)=beta_1
         l2=l2-1
1     continue
c     write(*,'(4f10.5)') (alpha(jot),jot=left-k+1,left)
      return
      end

