SUBROUTINE bcucof(y,y1,y2,y12,d1,d2,c)
	
	IMPLICIT NONE
	REAL, INTENT(IN) :: d1,d2
	REAL, DIMENSION(4), INTENT(IN) :: y,y1,y2,y12
	REAL, DIMENSION(4,4), INTENT(OUT) :: c
	REAL, DIMENSION(16) :: x
	REAL, DIMENSION(16,16) :: wt
	DATA wt /1,0,-3,2,4*0,-3,0,9,-6,2,0,-6,4,&
		8*0,3,0,-9,6,-2,0,6,-4,10*0,9,-6,2*0,-6,4,2*0,3,-2,6*0,-9,6,&
		2*0,6,-4,4*0,1,0,-3,2,-2,0,6,-4,1,0,-3,2,8*0,-1,0,3,-2,1,0,-3,&
		2,10*0,-3,2,2*0,3,-2,6*0,3,-2,2*0,-6,4,2*0,3,-2,0,1,-2,1,5*0,&
		-3,6,-3,0,2,-4,2,9*0,3,-6,3,0,-2,4,-2,10*0,-3,3,2*0,2,-2,2*0,&
		-1,1,6*0,3,-3,2*0,-2,2,5*0,1,-2,1,0,-2,4,-2,0,1,-2,1,9*0,-1,2,&
		-1,0,1,-2,1,10*0,1,-1,2*0,-1,1,6*0,-1,1,2*0,2,-2,2*0,-1,1/
	x(1:4)=y
	x(5:8)=y1*d1
	x(9:12)=y2*d2
	x(13:16)=y12*d1*d2
	x=matmul(wt,x)
	c=reshape(x,(/4,4/),order=(/2,1/))
	END SUBROUTINE bcucof