c
c
c***********************************************************************
c
      subroutine ORFAC2 ( kgp1, mmax, nu, b )
c
c-----------------------------------------------------------------------
c
c     stage 2:  special householder reduction of
c
c         kg     ( dr' . r3 )      (dr'' . r5 )
c                (-----. -- )      (-----. -- )
c     m-n-kg     (  0  . r4 )  to  (  0  . r4 )
c                (-----. -- )      (-----. -- )
c         kg     (nu*d . 0  )      (  0  . r6 )
c
c                   kg    1          kg    1
c
c     where dr', r3, and r4 are as in orfac1, nu is the
c     marquardt parameter, d is a diagonal matrix 
c     consisting of the lengths of the columns of dr',
c     and dr'' is in upper triangular form.  details in
c     (1), pp. 423-424.  note that the (m-n-kg) band of
c     zeroes, and r4, are omitted in storage.
c
c  ......................................................
c
      real b(mmax,kgp1),nu
c
      kg=kgp1-1
      kg2=2*kg
      kg23=kg2+3
      do 1070 k=1,kg
         kp1=k+1
         kgpk=kg+k
         kgpkm1=kgpk-1
         b(kgpk,k)=nu*b(kg23,k)
         alpha=b(k,k)**2
         do 1050 i=kgp1,kgpk
1050        alpha=alpha+b(i,k)**2
         alpha=sign(sqrt(alpha),b(k,k)) 
         u1=b(k,k)+alpha
         beta=alpha*u1
         b(k,k)=-alpha
c              the k-th reflection modifies only rows k,
c              kg+1, kg+2, ..., kg+k, and columns k to kg+1.
         do 1070 j=kp1,kgp1
            b(kgpk,j)=0.
            acum=u1*b(k,j)
            do 1060 i=kgp1,kgpkm1
1060           acum=acum+b(i,k)*b(i,j)
            acum=acum/beta
            b(k,j)=b(k,j)-u1*acum
            do 1070 i=kgp1,kgpk
1070           b(i,j)=b(i,j)-b(i,k)*acum
c
      return
c
      end 
