subroutine phi_tremble(residue,numberofatoms,dtheta,CAatomnumber,coord)
! finds phi1,psi1,phi2,psi2 with RRRR ~ 1
use defs; use noses; implicit none
integer(i4b), intent(in) :: residue,numberofatoms
real(dp), intent(in) :: dtheta  ! actual random angle scale
integer(i4b), dimension(:), intent(in) :: CAatomnumber
real(dp), dimension(:,:), intent(inout) :: coord 
! internal variables:
integer(i4b) :: amino
real(dp),dimension(3) :: N1,CA1,C1,N2,CA2,C2,axis1,axis2, &
                         axis3,axis4
real(dp),dimension(4) :: angles
real(dp),dimension(3,3) :: Axes
real(dp) :: norm,maxfactor,dphi,dpsi

! STUFF FOR DGESV
integer(i4b), parameter :: N = 3,NHRS=1,LDA=3,LDB=3
integer(i4b),dimension(3) :: IPIV
real(dp),dimension(3,1) :: B
integer(i4b) :: INFO 

N1 = coord(:,CAatomnumber(residue)-1)
CA1 = coord(:,CAatomnumber(residue))
C1 = coord(:,CAatomnumber(residue)+1) 
N2 = coord(:,CAatomnumber(residue+1)-1) 
CA2 = coord(:,CAatomnumber(residue+1))
C2 = coord(:,CAatomnumber(residue+1)+1)

! First find the values of the four small angles
axis1 = CA1 - N1; axis1 = axis1/sqrt(dot_product(axis1,axis1)) 
axis2 = C1 - CA1; axis2 = axis2/sqrt(dot_product(axis2,axis2))
axis3 = CA2 - N2; axis3 = axis3/sqrt(dot_product(axis3,axis3))
axis4 = C2 - CA2; axis4 = axis4/sqrt(dot_product(axis4,axis4))  
Axes(:,1) = axis1; Axes(:,2) = axis2;  Axes(:,3) = axis3
B(:,1) = axis4 
! DGESV finds X such that Axes * X = B
call DGESV ( N, NHRS, Axes, LDA, IPIV, B, LDB, INFO )
! on exit, B = X.
angles(1) = B(1,1); angles(2) = B(2,1)
angles(3) = B(3,1); angles(4) = -1.0_dp
norm = dot_product(angles,angles)
maxfactor = dtheta/sqrt(norm)
angles = maxfactor*angles
dphi = angles(1); amino = residue  
call phi_rot(amino,coord,dphi,Caatomnumber,numberofatoms)
dpsi = angles(2); amino = residue
call psi_rot(amino,coord,dpsi,Caatomnumber,numberofatoms)
dphi = angles(3); amino = residue + 1
call phi_rot(amino,coord,dphi,Caatomnumber,numberofatoms)
dpsi = angles(4); amino = residue + 1
call psi_rot(amino,coord,dpsi,Caatomnumber,numberofatoms)

end subroutine phi_tremble