1 C Output from Public domain Ratfor, version 1.0
2 subroutine stxwx(x,z,w,k,xknot,n,y,hs0,hs1,hs2,hs3)
6 DOUBLE precision x(k),z(k),w(k), xknot(n+4),y(n),
7 & hs0(n),hs1(n),hs2(n),hs3(n)
9 DOUBLE precision eps,vnikx(4,1),work(16)
10 integer lenxk, i,j, ileft,mflag
16 C Initialise the output vectors
25 C Compute X' W^2 X -> hs0,hs1,hs2,hs3 and X' W^2 Z -> y
26 C Note that here the weights w(i) == sqrt(wt[i]) where wt[] where original weights
31 call interv(xknot(1), n+1, x(i), ileft, mflag)
32 C if(mflag==-1) {write(6,'("Error in hess ",i2)')mflag;stop}
33 C if(mflag==-1) {return}
35 if(x(i).le.(xknot(ileft)+eps))then
40 C else{write(6,'("Error in hess ",i2)')mflag;stop}}
42 call bsplvd (xknot,lenxk,4,x(i),ileft,work,vnikx,1)
45 y(j) = y(j)+w(i)**2*z(i)*vnikx(1,1)
46 hs0(j)=hs0(j)+w(i)**2*vnikx(1,1)**2
47 hs1(j)=hs1(j)+w(i)**2*vnikx(1,1)*vnikx(2,1)
48 hs2(j)=hs2(j)+w(i)**2*vnikx(1,1)*vnikx(3,1)
49 hs3(j)=hs3(j)+w(i)**2*vnikx(1,1)*vnikx(4,1)
51 y(j) = y(j)+w(i)**2*z(i)*vnikx(2,1)
52 hs0(j)=hs0(j)+w(i)**2*vnikx(2,1)**2
53 hs1(j)=hs1(j)+w(i)**2*vnikx(2,1)*vnikx(3,1)
54 hs2(j)=hs2(j)+w(i)**2*vnikx(2,1)*vnikx(4,1)
56 y(j) = y(j)+w(i)**2*z(i)*vnikx(3,1)
57 hs0(j)=hs0(j)+w(i)**2*vnikx(3,1)**2
58 hs1(j)=hs1(j)+w(i)**2*vnikx(3,1)*vnikx(4,1)
60 y(j) = y(j)+w(i)**2*z(i)*vnikx(4,1)
61 hs0(j)=hs0(j)+w(i)**2*vnikx(4,1)**2