subroutine stitch(hbp,hbt,hba,hrp,hrt,hra,dp2,dt2)
*
*  stitch together the envelope and the interior. 
*
      implicit double precision(a-h,o-z)

      common/tablec/tablec(6,39,87)
      common/tableo/tableo(6,36,87)
      common/shells/shells(8000)
      common/temp/s1,r1,s2,r2,b2,p2,t2,e2,xc2,xo2,fca2,f2,q2,w2,c
      common/surf/u(2,3),v(2,3),w(2,3),r,b,i,k,l,m,kstart,npass

      dimension w1(2,2),ww2(2,2),w3(2,2),w4(2,2),h(2,2),del(2),kk(3)
      data fnat/0.4342945/
*
* initialize variables
*
      b=dlog10(b2)
      do 1 i=1,3
        kk(i)=0
1     continue
*
*  start loop where we try to trap luminosity and radius.
*  2 is the top of the loop.
*
2     g1=k*((w(1,2)-w(1,3))*(r-w(2,2))+(w(2,3)-w(2,2))*(b-w(1,2)))
      if ( g1 .lt. 0. ) then
        do 4 i=1,2 
          w(i,1)=w(i,3)+w(i,2)-w(i,1) 
4       continue
        k=-k
        m=m-3*l
        l=-l
        kk(1)=1
        go to 2
      endif
      g2=k*((w(1,3)-w(1,1))*(r-w(2,3))+(w(2,1)-w(2,3))*(b-w(1,3)))
      if ( g2 .lt. 0. ) then
         do 7 i=1,2 
            w(i,2)=2*w(i,1)-w(i,2)
7        continue
         k=-k
         kk(2)=1
         go to 2
      endif
      g3=k*((w(1,1)-w(1,2))*(r-w(2,1))+(w(2,2)-w(2,1))*(b-w(1,1)))
      if ( g3 .lt. 0. ) then
         do 10 i=1,2
            w(i,3)=2*w(i,1)-w(i,3)
10       continue
         k=-k
         kk(3)=1
         go to 2
      endif
*
* this is the bottom of the loop
* now that we've trapped the values, let's get on with it
*
      if ( kstart .eq. 1 ) then
         kk(1)=1
         kk(2)=1
         kk(3)=1
         kstart=0
      endif
      if ( kk(1)+kk(2)+kk(3) .gt. 0 ) then
         write(11) shells
         rewind 11
         do 13 i=1,3
            if (kk(i) .gt. 0) call env
13       continue
         read(11) shells
         rewind 11
      elseif ( npass .gt. 0 ) then
         go to 17
      endif
*
* goes through here just after starting up the sequence. 
*
      do 15 i=1,2
         do 25 j=1,2
            n=m-l*j 
            cd=w(j,n)-w(j,1)
            w1(i,j)=(u(i,n)-u(i,1))/cd
            ww2(i,j)=(v(i,n)-v(i,1))/cd 
25       continue
15    continue

      call solve(w3(1,1),w3(2,1),w1,-1.0d0,0.0d0)
      call solve(w3(1,2),w3(2,2),w1,0.0d0,-1.0d0)

      do 16 i=1,2
         do 26 j=1,2
            w4(i,j)=ww2(i,1)*w3(1,j)+ww2(i,2)*w3(2,j)
26       continue
16    continue
17    del(1)=p2-u(1,1)
      del(2)=t2-u(2,1)
      r=w(2,1)
      h1=v(1,1)
      h2=v(2,1)
      do 18 i=1,2
         r=r+w3(2,i)*del(i)
         h1=h1+w4(1,i)*del(i) 
         h2=h2+w4(2,i)*del(i) 
18    continue
      h1=hba+b2-10.**h1
      h2=hra+r2-h2
      h(1,1)=hbp-w4(1,1)*b2/fnat
      h(1,2)=hbt-w4(1,2)*b2/fnat
      h(2,1)=hrp-w4(2,1)
      h(2,2)=hrt-w4(2,2)
      call solve(dp2,dt2,h,h1,h2)
      npass=npass+1 

      return
      end 

************************************************************************