src_CSA_DiL removed from prerelease, current version in devel
[unres.git] / source / unres / src_CSA_DiL / local_move.f
diff --git a/source/unres/src_CSA_DiL/local_move.f b/source/unres/src_CSA_DiL/local_move.f
deleted file mode 100644 (file)
index 763d3cc..0000000
+++ /dev/null
@@ -1,972 +0,0 @@
-c-------------------------------------------------------------
-
-      subroutine local_move_init(debug)
-crc      implicit none
-
-c     Includes
-      implicit real*8 (a-h,o-z)
-      include 'DIMENSIONS'  ! Needed by COMMON.LOCAL
-      include 'COMMON.GEO'  ! For pi, deg2rad
-      include 'COMMON.LOCAL'  ! For vbl
-      include 'COMMON.LOCMOVE'
-
-c     INPUT arguments
-      logical debug
-
-
-c     Determine wheter to do some debugging output
-      locmove_output=debug
-
-c     Set the init_called flag to 1
-      init_called=1
-
-c     The following are never changed
-      min_theta=60.D0*deg2rad  ! (0,PI)
-      max_theta=175.D0*deg2rad  ! (0,PI)
-      dmin2=vbl*vbl*2.*(1.-cos(min_theta))
-      dmax2=vbl*vbl*2.*(1.-cos(max_theta))
-      flag=1.0D300
-      small=1.0D-5
-      small2=0.5*small*small
-
-c     Not really necessary...
-      a_n=0
-      b_n=0
-      res_n=0
-
-      return
-      end
-
-c-------------------------------------------------------------
-
-      subroutine local_move(n_start, n_end, PHImin, PHImax)
-c     Perform a local move between residues m and n (inclusive)
-c     PHImin and PHImax [0,PI] determine the size of the move
-c     Works on whatever structure is in the variables theta and phi,
-c     sidechain variables are left untouched
-c     The final structure is NOT minimized, but both the cartesian
-c     variables c and the angles are up-to-date at the end (no further
-c     chainbuild is required)
-crc      implicit none
-
-c     Includes
-      implicit real*8 (a-h,o-z)
-      include 'DIMENSIONS'
-      include 'COMMON.GEO'
-      include 'COMMON.CHAIN'
-      include 'COMMON.VAR'
-      include 'COMMON.MINIM'
-      include 'COMMON.SBRIDGE'
-      include 'COMMON.LOCMOVE'
-
-c     External functions
-      integer move_res
-      external move_res
-      double precision ran_number
-      external ran_number
-
-c     INPUT arguments
-      integer n_start, n_end  ! First and last residues to move
-      double precision PHImin, PHImax  ! min/max angles [0,PI]
-
-c     Local variables
-      integer i,j
-      double precision min,max
-      integer iretcode
-
-
-c     Check if local_move_init was called.  This assumes that it
-c     would not be 1 if not explicitely initialized
-      if (init_called.ne.1) then
-        write(6,*)'   ***   local_move_init not called!!!'
-        stop
-      endif
-
-c     Quick check for crazy range
-      if (n_start.gt.n_end .or. n_start.lt.1 .or. n_end.gt.nres) then
-        write(6,'(a,i3,a,i3)')
-     +       '   ***   Cannot make local move between n_start = ',
-     +       n_start,' and n_end = ',n_end
-        return
-      endif
-
-c     Take care of end residues first...
-      if (n_start.eq.1) then
-c     Move residue 1 (completely random)
-        theta(3)=ran_number(min_theta,max_theta)
-        phi(4)=ran_number(-PI,PI)
-        i=2
-      else
-        i=n_start
-      endif
-      if (n_end.eq.nres) then
-c     Move residue nres (completely random)
-        theta(nres)=ran_number(min_theta,max_theta)
-        phi(nres)=ran_number(-PI,PI)
-        j=nres-1
-      else
-        j=n_end
-      endif
-
-c     ...then go through all other residues one by one
-c     Start from the two extremes and converge
-      call chainbuild
-      do while (i.le.j)
-        min=PHImin
-        max=PHImax
-c$$$c     Move the first two residues by less than the others
-c$$$        if (i-n_start.lt.3) then
-c$$$          if (i-n_start.eq.0) then
-c$$$            min=0.4*PHImin
-c$$$            max=0.4*PHImax
-c$$$          else if (i-n_start.eq.1) then
-c$$$            min=0.8*PHImin
-c$$$            max=0.8*PHImax
-c$$$          else if (i-n_start.eq.2) then
-c$$$            min=PHImin
-c$$$            max=PHImax
-c$$$          endif
-c$$$        endif
-
-c     The actual move, on residue i
-        iretcode=move_res(min,max,i)  ! Discard iretcode
-        i=i+1
-
-        if (i.le.j) then
-          min=PHImin
-          max=PHImax
-c$$$c     Move the last two residues by less than the others
-c$$$          if (n_end-j.lt.3) then
-c$$$            if (n_end-j.eq.0) then
-c$$$              min=0.4*PHImin
-c$$$              max=0.4*PHImax
-c$$$            else if (n_end-j.eq.1) then
-c$$$              min=0.8*PHImin
-c$$$              max=0.8*PHImax
-c$$$            else if (n_end-j.eq.2) then
-c$$$              min=PHImin
-c$$$              max=PHImax
-c$$$            endif
-c$$$          endif
-
-c     The actual move, on residue j
-          iretcode=move_res(min,max,j)  ! Discard iretcode
-          j=j-1
-        endif
-      enddo
-
-      call int_from_cart(.false.,.false.)
-
-      return
-      end
-
-c-------------------------------------------------------------
-
-      subroutine output_tabs
-c     Prints out the contents of a_..., b_..., res_...
-      implicit none
-
-c     Includes
-      include 'COMMON.GEO'
-      include 'COMMON.LOCMOVE'
-
-c     Local variables
-      integer i,j
-
-
-      write(6,*)'a_...'
-      write(6,'(8f7.1)')(a_ang(i)*rad2deg,i=0,a_n-1)
-      write(6,'(8(2x,3l1,2x))')((a_tab(i,j),i=0,2),j=0,a_n-1)
-
-      write(6,*)'b_...'
-      write(6,'(4f7.1)')(b_ang(i)*rad2deg,i=0,b_n-1)
-      write(6,'(4(2x,3l1,2x))')((b_tab(i,j),i=0,2),j=0,b_n-1)
-
-      write(6,*)'res_...'
-      write(6,'(12f7.1)')(res_ang(i)*rad2deg,i=0,res_n-1)
-      write(6,'(12(2x,3l1,2x))')((res_tab(0,i,j),i=0,2),j=0,res_n-1)
-      write(6,'(12(2x,3l1,2x))')((res_tab(1,i,j),i=0,2),j=0,res_n-1)
-      write(6,'(12(2x,3l1,2x))')((res_tab(2,i,j),i=0,2),j=0,res_n-1)
-
-      return
-      end
-
-c-------------------------------------------------------------
-
-      subroutine angles2tab(PHImin,PHImax,n,ang,tab)
-c     Only uses angles if [0,PI] (but PHImin cannot be 0.,
-c     and PHImax cannot be PI)
-      implicit none
-
-c     Includes
-      include 'COMMON.GEO'
-
-c     INPUT arguments
-      double precision PHImin,PHImax
-
-c     OUTPUT arguments
-      integer n
-      double precision ang(0:3)
-      logical tab(0:2,0:3)
-
-
-      if (PHImin .eq. PHImax) then
-c     Special case with two 010's
-        n = 2;
-        ang(0) = -PHImin;
-        ang(1) = PHImin;
-        tab(0,0) = .false.
-        tab(2,0) = .false.
-        tab(0,1) = .false.
-        tab(2,1) = .false.
-        tab(1,0) = .true.
-        tab(1,1) = .true.
-      else if (PHImin .eq. PI) then
-c     Special case with one 010
-        n = 1
-        ang(0) = PI
-        tab(0,0) = .false.
-        tab(2,0) = .false.
-        tab(1,0) = .true.
-      else if (PHImax .eq. 0.) then
-c     Special case with one 010
-        n = 1
-        ang(0) = 0.
-        tab(0,0) = .false.
-        tab(2,0) = .false.
-        tab(1,0) = .true.
-      else
-c     Standard cases
-        n = 0
-        if (PHImin .gt. 0.) then
-c     Start of range (011)
-          ang(n) = PHImin
-          tab(0,n) = .false.
-          tab(1,n) = .true.
-          tab(2,n) = .true.
-c     End of range (110)
-          ang(n+1) = -PHImin
-          tab(0,n+1) = .true.
-          tab(1,n+1) = .true.
-          tab(2,n+1) = .false.
-          n = n+2
-        endif
-        if (PHImax .lt. PI) then
-c     Start of range (011)
-          ang(n) = -PHImax
-          tab(0,n) = .false.
-          tab(1,n) = .true.
-          tab(2,n) = .true.
-c     End of range (110)
-          ang(n+1) = PHImax
-          tab(0,n+1) = .true.
-          tab(1,n+1) = .true.
-          tab(2,n+1) = .false.
-          n = n+2
-        endif
-      endif
-
-      return
-      end
-
-c-------------------------------------------------------------
-
-      subroutine minmax_angles(x,y,z,r,n,ang,tab)
-c     When solutions do not exist, assume all angles
-c     are acceptable - i.e., initial geometry must be correct
-      implicit none
-
-c     Includes
-      include 'COMMON.GEO'
-      include 'COMMON.LOCMOVE'
-
-c     Input arguments
-      double precision x,y,z,r
-
-c     Output arguments
-      integer n
-      double precision ang(0:3)
-      logical tab(0:2,0:3)
-
-c     Local variables
-      double precision num, denom, phi
-      double precision Kmin, Kmax
-      integer i
-
-
-      num = x*x + y*y + z*z
-      denom = x*x + y*y
-      n = 0
-      if (denom .gt. 0.) then
-        phi = atan2(y,x)
-        denom = 2.*r*sqrt(denom)
-        num = num+r*r
-        Kmin = (num - dmin2)/denom
-        Kmax = (num - dmax2)/denom
-
-c     Allowed values of K (else all angles are acceptable)
-c     -1 <= Kmin <  1
-c     -1 <  Kmax <= 1
-        if (Kmin .gt. 1. .or. abs(Kmin-1.) .lt. small2) then
-          Kmin = -flag
-        else if (Kmin .lt. -1. .or. abs(Kmin+1.) .lt. small2) then
-          Kmin = PI
-        else
-          Kmin = acos(Kmin)
-        endif
-
-        if (Kmax .lt. -1. .or. abs(Kmax+1.) .lt. small2) then
-          Kmax = flag
-        else if (Kmax .gt. 1. .or. abs(Kmax-1.) .lt. small2) then
-          Kmax = 0.
-        else
-          Kmax = acos(Kmax)
-        endif
-
-        if (Kmax .lt. Kmin) Kmax = Kmin
-
-        call angles2tab(Kmin, Kmax, n, ang, tab)
-
-c     Add phi and check that angles are within range (-PI,PI]
-        do i=0,n-1
-          ang(i) = ang(i)+phi
-          if (ang(i) .le. -PI) then
-            ang(i) = ang(i)+2.*PI
-          else if (ang(i) .gt. PI) then
-            ang(i) = ang(i)-2.*PI
-          endif
-        enddo
-      endif
-
-      return
-      end
-
-c-------------------------------------------------------------
-
-      subroutine construct_tab
-c     Take a_... and b_... values and produces the results res_...
-c     x_ang are assumed to be all different (diff > small)
-c     x_tab(1,i) must be 1 for all i (i.e., all x_ang are acceptable)
-      implicit none
-
-c     Includes
-      include 'COMMON.LOCMOVE'
-
-c     Local variables
-      integer n_max,i,j,index
-      logical done
-      double precision phi
-
-
-      n_max = a_n + b_n
-      if (n_max .eq. 0) then
-        res_n = 0
-        return
-      endif
-
-      do i=0,n_max-1
-        do j=0,1
-          res_tab(j,0,i) = .true.
-          res_tab(j,2,i) = .true.
-          res_tab(j,1,i) = .false.
-        enddo
-      enddo
-
-      index = 0
-      phi = -flag
-      done = .false.
-      do while (.not.done)
-        res_ang(index) = flag
-
-c     Check a first...
-        do i=0,a_n-1
-          if ((a_ang(i)-phi).gt.small .and.
-     +         a_ang(i) .lt. res_ang(index)) then
-c     Found a lower angle
-            res_ang(index) = a_ang(i)
-c     Copy the values from a_tab into res_tab(0,,)
-            res_tab(0,0,index) = a_tab(0,i)
-            res_tab(0,1,index) = a_tab(1,i)
-            res_tab(0,2,index) = a_tab(2,i)
-c     Set default values for res_tab(1,,)
-            res_tab(1,0,index) = .true.
-            res_tab(1,1,index) = .false.
-            res_tab(1,2,index) = .true.
-          else if (abs(a_ang(i)-res_ang(index)).lt.small) then
-c     Found an equal angle (can only be equal to a b_ang)
-            res_tab(0,0,index) = a_tab(0,i)
-            res_tab(0,1,index) = a_tab(1,i)
-            res_tab(0,2,index) = a_tab(2,i)
-          endif
-        enddo
-c     ...then check b
-        do i=0,b_n-1
-          if ((b_ang(i)-phi).gt.small .and.
-     +         b_ang(i) .lt. res_ang(index)) then
-c     Found a lower angle
-            res_ang(index) = b_ang(i)
-c     Copy the values from b_tab into res_tab(1,,)
-            res_tab(1,0,index) = b_tab(0,i)
-            res_tab(1,1,index) = b_tab(1,i)
-            res_tab(1,2,index) = b_tab(2,i)
-c     Set default values for res_tab(0,,)
-            res_tab(0,0,index) = .true.
-            res_tab(0,1,index) = .false.
-            res_tab(0,2,index) = .true.
-          else if (abs(b_ang(i)-res_ang(index)).lt.small) then
-c     Found an equal angle (can only be equal to an a_ang)
-            res_tab(1,0,index) = b_tab(0,i)
-            res_tab(1,1,index) = b_tab(1,i)
-            res_tab(1,2,index) = b_tab(2,i)
-          endif
-        enddo
-
-        if (res_ang(index) .eq. flag) then
-          res_n = index
-          done = .true.
-        else if (index .eq. n_max-1) then
-          res_n = n_max
-          done = .true.
-        else
-          phi = res_ang(index)  ! Store previous angle
-          index = index+1
-        endif
-      enddo
-
-c     Fill the gaps
-c     First a...
-      index = 0
-      if (a_n .gt. 0) then
-        do while (.not.res_tab(0,1,index))
-          index=index+1
-        enddo
-        done = res_tab(0,2,index)
-        do i=index+1,res_n-1
-          if (res_tab(0,1,i)) then
-            done = res_tab(0,2,i)
-          else
-            res_tab(0,0,i) = done
-            res_tab(0,1,i) = done
-            res_tab(0,2,i) = done
-          endif
-        enddo
-        done = res_tab(0,0,index)
-        do i=index-1,0,-1
-          if (res_tab(0,1,i)) then
-            done = res_tab(0,0,i)
-          else
-            res_tab(0,0,i) = done
-            res_tab(0,1,i) = done
-            res_tab(0,2,i) = done
-          endif
-        enddo
-      else
-        do i=0,res_n-1
-          res_tab(0,0,i) = .true.
-          res_tab(0,1,i) = .true.
-          res_tab(0,2,i) = .true.
-        enddo
-      endif
-c     ...then b
-      index = 0
-      if (b_n .gt. 0) then
-        do while (.not.res_tab(1,1,index))
-          index=index+1
-        enddo
-        done = res_tab(1,2,index)
-        do i=index+1,res_n-1
-          if (res_tab(1,1,i)) then
-            done = res_tab(1,2,i)
-          else
-            res_tab(1,0,i) = done
-            res_tab(1,1,i) = done
-            res_tab(1,2,i) = done
-          endif
-        enddo
-        done = res_tab(1,0,index)
-        do i=index-1,0,-1
-          if (res_tab(1,1,i)) then
-            done = res_tab(1,0,i)
-          else
-            res_tab(1,0,i) = done
-            res_tab(1,1,i) = done
-            res_tab(1,2,i) = done
-          endif
-        enddo
-      else
-        do i=0,res_n-1
-          res_tab(1,0,i) = .true.
-          res_tab(1,1,i) = .true.
-          res_tab(1,2,i) = .true.
-        enddo
-      endif
-
-c     Finally fill the last row with AND operation
-      do i=0,res_n-1
-        do j=0,2
-          res_tab(2,j,i) = (res_tab(0,j,i) .and. res_tab(1,j,i))
-        enddo
-      enddo
-
-      return 
-      end
-
-c-------------------------------------------------------------
-
-      subroutine construct_ranges(phi_n,phi_start,phi_end)
-c     Given the data in res_..., construct a table of 
-c     min/max allowed angles
-      implicit none
-
-c     Includes
-      include 'COMMON.GEO'
-      include 'COMMON.LOCMOVE'
-
-c     Output arguments
-      integer phi_n
-      double precision phi_start(0:11),phi_end(0:11)
-
-c     Local variables
-      logical done
-      integer index
-
-
-      if (res_n .eq. 0) then
-c     Any move is allowed
-        phi_n = 1
-        phi_start(0) = -PI
-        phi_end(0) = PI
-      else
-        phi_n = 0
-        index = 0
-        done = .false.
-        do while (.not.done)
-c     Find start of range (01x)
-          done = .false.
-          do while (.not.done)
-            if (res_tab(2,0,index).or.(.not.res_tab(2,1,index))) then
-              index=index+1
-            else
-              done = .true.
-              phi_start(phi_n) = res_ang(index)
-            endif
-            if (index .eq. res_n) done = .true.
-          enddo
-c     If a start was found (index < res_n), find the end of range (x10)
-c     It may not be found without wrapping around
-          if (index .lt. res_n) then
-            done = .false.
-            do while (.not.done)
-              if ((.not.res_tab(2,1,index)).or.res_tab(2,2,index)) then
-                index=index+1
-              else
-                done = .true.
-              endif
-              if (index .eq. res_n) done = .true.
-            enddo
-            if (index .lt. res_n) then
-c     Found the end of the range
-              phi_end(phi_n) = res_ang(index)
-              phi_n=phi_n+1
-              index=index+1
-              if (index .eq. res_n) then
-                done = .true.
-              else
-                done = .false.
-              endif
-            else
-c     Need to wrap around
-              done = .true.
-              phi_end(phi_n) = flag
-            endif
-          endif
-        enddo
-c     Take care of the last one if need to wrap around
-        if (phi_end(phi_n) .eq. flag) then
-          index = 0
-          do while ((.not.res_tab(2,1,index)).or.res_tab(2,2,index))
-            index=index+1
-          enddo
-          phi_end(phi_n) = res_ang(index) + 2.*PI
-          phi_n=phi_n+1
-        endif
-      endif
-
-      return
-      end
-
-c-------------------------------------------------------------
-
-      subroutine fix_no_moves(phi)
-      implicit none
-
-c     Includes
-      include 'COMMON.GEO'
-      include 'COMMON.LOCMOVE'
-
-c     Output arguments
-      double precision phi
-
-c     Local variables
-      integer index
-      double precision diff,temp
-
-
-c     Look for first 01x in gammas (there MUST be at least one)
-      diff = flag
-      index = 0
-      do while (res_tab(1,0,index) .or. (.not.res_tab(1,1,index)))
-        index=index+1
-      enddo
-      if (res_ang(index) .le. 0.D0) then ! Make sure it's from PHImax
-c     Try to increase PHImax
-        if (index .gt. 0) then
-          phi = res_ang(index-1)
-          diff = abs(res_ang(index) - res_ang(index-1))
-        endif
-c     Look for last (corresponding) x10
-        index = res_n - 1
-        do while ((.not.res_tab(1,1,index)) .or. res_tab(1,2,index))
-          index=index-1
-        enddo
-        if (index .lt. res_n-1) then
-          temp = abs(res_ang(index) - res_ang(index+1))
-          if (temp .lt. diff) then
-            phi = res_ang(index+1)
-            diff = temp
-          endif
-        endif
-      endif
-
-c     If increasing PHImax didn't work, decreasing PHImin
-c     will (with one exception)
-c     Look for first x10 (there MUST be at least one)
-      index = 0
-      do while ((.not.res_tab(1,1,index)) .or. res_tab(1,2,index))
-        index=index+1
-      enddo
-      if (res_ang(index) .lt. 0.D0) then ! Make sure it's from PHImin
-c     Try to decrease PHImin
-        if (index .lt. res_n-1) then
-          temp = abs(res_ang(index) - res_ang(index+1))
-          if (res_ang(index+1) .le. 0.D0 .and. temp .lt. diff) then
-            phi = res_ang(index+1)
-            diff = temp
-          endif
-        endif
-c     Look for last (corresponding) 01x
-        index = res_n - 1
-        do while (res_tab(1,0,index) .or. (.not.res_tab(1,1,index)))
-          index=index-1
-        enddo
-        if (index .gt. 0) then
-          temp = abs(res_ang(index) - res_ang(index-1))
-          if (res_ang(index-1) .ge. 0.D0 .and. temp .lt. diff) then
-            phi = res_ang(index-1)
-            diff = temp
-          endif
-        endif
-      endif
-
-c     If it still didn't work, it must be PHImax == 0. or PHImin == PI
-      if (diff .eq. flag) then
-        index = 0
-        if (res_tab(index,1,0) .or. (.not.res_tab(index,1,1)) .or.
-     +       res_tab(index,1,2)) index = res_n - 1
-c     This MUST work at this point
-        if (index .eq. 0) then
-          phi = res_ang(1)
-        else
-          phi = res_ang(index - 1)
-        endif
-      endif
-
-      return
-      end
-
-c-------------------------------------------------------------
-
-      integer function move_res(PHImin,PHImax,i_move)
-c     Moves residue i_move (in array c), leaving everything else fixed
-c     Starting geometry is not checked, it should be correct!
-c     R(,i_move) is the only residue that will move, but must have
-c     1 < i_move < nres (i.e., cannot move ends)
-c     Whether any output is done is controlled by locmove_output
-crc      implicit none
-
-c     Includes
-      implicit real*8 (a-h,o-z)
-      include 'DIMENSIONS'
-      include 'COMMON.CHAIN'
-      include 'COMMON.GEO'
-      include 'COMMON.LOCMOVE'
-
-c     External functions
-      double precision ran_number
-      external ran_number
-
-c     Input arguments
-      double precision PHImin,PHImax
-      integer i_move
-
-c     RETURN VALUES:
-c     0: move successfull
-c     1: Dmin or Dmax had to be modified
-c     2: move failed - check your input geometry
-
-
-c     Local variables
-      double precision X(0:2),Y(0:2),Z(0:2),Orig(0:2)
-      double precision P(0:2)
-      logical no_moves,done
-      integer index,i,j
-      double precision phi,temp,radius
-      double precision phi_start(0:11), phi_end(0:11)
-      integer phi_n
-
-c     Set up the coordinate system
-      do i=0,2
-        Orig(i)=0.5*(c(i+1,i_move-1)+c(i+1,i_move+1)) ! Position of origin
-      enddo
-
-      do i=0,2
-        Z(i)=c(i+1,i_move+1)-c(i+1,i_move-1)
-      enddo
-      temp=sqrt(Z(0)*Z(0)+Z(1)*Z(1)+Z(2)*Z(2))
-      do i=0,2
-        Z(i)=Z(i)/temp
-      enddo
-
-      do i=0,2
-        X(i)=c(i+1,i_move)-Orig(i)
-      enddo
-c     radius is the radius of the circle on which c(,i_move) can move
-      radius=sqrt(X(0)*X(0)+X(1)*X(1)+X(2)*X(2))
-      do i=0,2
-        X(i)=X(i)/radius
-      enddo
-
-      Y(0)=Z(1)*X(2)-X(1)*Z(2)
-      Y(1)=X(0)*Z(2)-Z(0)*X(2)
-      Y(2)=Z(0)*X(1)-X(0)*Z(1)
-
-c     Calculate min, max angles coming from dmin, dmax to c(,i_move-2)
-      if (i_move.gt.2) then
-        do i=0,2
-          P(i)=c(i+1,i_move-2)-Orig(i)
-        enddo
-        call minmax_angles(P(0)*X(0)+P(1)*X(1)+P(2)*X(2),
-     +       P(0)*Y(0)+P(1)*Y(1)+P(2)*Y(2),
-     +       P(0)*Z(0)+P(1)*Z(1)+P(2)*Z(2),
-     +       radius,a_n,a_ang,a_tab)
-      else
-        a_n=0
-      endif
-
-c     Calculate min, max angles coming from dmin, dmax to c(,i_move+2)
-      if (i_move.lt.nres-2) then
-        do i=0,2
-          P(i)=c(i+1,i_move+2)-Orig(i)
-        enddo
-        call minmax_angles(P(0)*X(0)+P(1)*X(1)+P(2)*X(2),
-     +       P(0)*Y(0)+P(1)*Y(1)+P(2)*Y(2),
-     +       P(0)*Z(0)+P(1)*Z(1)+P(2)*Z(2),
-     +       radius,b_n,b_ang,b_tab)
-      else
-        b_n=0
-      endif
-
-c     Construct the resulting table for alpha and beta
-      call construct_tab()
-
-      if (locmove_output) then
-        print *,'ALPHAS & BETAS TABLE'
-        call output_tabs()
-      endif
-
-c     Check that there is at least one possible move
-      no_moves = .true.
-      if (res_n .eq. 0) then
-        no_moves = .false.
-      else
-        index = 0
-        do while ((index .lt. res_n) .and. no_moves)
-          if (res_tab(2,1,index)) no_moves = .false.
-          index=index+1
-        enddo
-      endif
-      if (no_moves) then
-        if (locmove_output) print *,'   ***   Cannot move anywhere'
-        move_res=2
-        return
-      endif
-
-c     Transfer res_... into a_...
-      a_n = 0
-      do i=0,res_n-1
-        if ( (res_tab(2,0,i).neqv.res_tab(2,1,i)) .or.
-     +       (res_tab(2,0,i).neqv.res_tab(2,2,i)) ) then
-          a_ang(a_n) = res_ang(i)
-          do j=0,2
-            a_tab(j,a_n) = res_tab(2,j,i)
-          enddo
-          a_n=a_n+1
-        endif
-      enddo
-
-c     Check that the PHI's are within [0,PI]
-      if (PHImin .lt. 0. .or. abs(PHImin) .lt. small) PHImin = -flag
-      if (PHImin .gt. PI .or. abs(PHImin-PI) .lt. small) PHImin = PI
-      if (PHImax .gt. PI .or. abs(PHImax-PI) .lt. small) PHImax = flag
-      if (PHImax .lt. 0. .or. abs(PHImax) .lt. small) PHImax = 0.
-      if (PHImax .lt. PHImin) PHImax = PHImin
-c     Calculate min and max angles coming from PHImin and PHImax,
-c     and put them in b_...
-      call angles2tab(PHImin, PHImax, b_n, b_ang, b_tab)
-c     Construct the final table
-      call construct_tab()
-
-      if (locmove_output) then
-        print *,'FINAL TABLE'
-        call output_tabs()
-      endif
-
-c     Check that there is at least one possible move
-      no_moves = .true.
-      if (res_n .eq. 0) then
-        no_moves = .false.
-      else
-        index = 0
-        do while ((index .lt. res_n) .and. no_moves)
-          if (res_tab(2,1,index)) no_moves = .false.
-          index=index+1
-        enddo
-      endif
-
-      if (no_moves) then
-c     Take care of the case where no solution exists...
-        call fix_no_moves(phi)
-        if (locmove_output) then
-          print *,'   ***   Had to modify PHImin or PHImax'
-          print *,'phi: ',phi*rad2deg
-        endif
-        move_res=1
-      else
-c     ...or calculate the solution
-c     Construct phi_start/phi_end arrays
-        call construct_ranges(phi_n, phi_start, phi_end)
-c     Choose random angle phi in allowed range(s)
-        temp = 0.
-        do i=0,phi_n-1
-          temp = temp + phi_end(i) - phi_start(i)
-        enddo
-        phi = ran_number(phi_start(0),phi_start(0)+temp)
-        index = 0
-        done = .false.
-        do while (.not.done)
-          if (phi .lt. phi_end(index)) then
-            done = .true.
-          else
-            index=index+1
-          endif
-          if (index .eq. phi_n) then
-            done = .true.
-          else if (.not.done) then
-            phi = phi + phi_start(index) - phi_end(index-1)
-          endif
-        enddo
-        if (index.eq.phi_n) phi=phi_end(phi_n-1) ! Fix numerical errors
-        if (phi .gt. PI) phi = phi-2.*PI
-
-        if (locmove_output) then
-          print *,'ALLOWED RANGE(S)'
-          do i=0,phi_n-1
-            print *,phi_start(i)*rad2deg,phi_end(i)*rad2deg
-          enddo
-          print *,'phi: ',phi*rad2deg
-        endif
-        move_res=0
-      endif
-
-c     Re-use radius as temp variable
-      temp=radius*cos(phi)
-      radius=radius*sin(phi)
-      do i=0,2
-        c(i+1,i_move)=Orig(i)+temp*X(i)+radius*Y(i)
-      enddo
-
-      return
-      end
-
-c-------------------------------------------------------------
-
-      subroutine loc_test
-crc      implicit none
-
-c     Includes
-      implicit real*8 (a-h,o-z)
-      include 'DIMENSIONS'
-      include 'COMMON.GEO'
-      include 'COMMON.LOCAL'
-      include 'COMMON.LOCMOVE'
-
-c     External functions
-      integer move_res
-      external move_res
-
-c     Local variables
-      integer i,j
-      integer phi_n
-      double precision phi_start(0:11),phi_end(0:11)
-      double precision phi
-      double precision R(0:2,0:5)
-
-      locmove_output=.true.
-
-c      call angles2tab(30.*deg2rad,70.*deg2rad,a_n,a_ang,a_tab)
-c      call angles2tab(80.*deg2rad,130.*deg2rad,b_n,b_ang,b_tab)
-c      call minmax_angles(0.D0,3.8D0,0.D0,3.8D0,b_n,b_ang,b_tab)
-c      call construct_tab
-c      call output_tabs
-
-c      call construct_ranges(phi_n,phi_start,phi_end)
-c      do i=0,phi_n-1
-c        print *,phi_start(i)*rad2deg,phi_end(i)*rad2deg
-c      enddo
-
-c      call fix_no_moves(phi)
-c      print *,'NO MOVES FOUND, BEST PHI IS',phi*rad2deg
-
-      R(0,0)=0.D0
-      R(1,0)=0.D0
-      R(2,0)=0.D0
-      R(0,1)=0.D0
-      R(1,1)=-cos(28.D0*deg2rad)
-      R(2,1)=-0.5D0-sin(28.D0*deg2rad)
-      R(0,2)=0.D0
-      R(1,2)=0.D0
-      R(2,2)=-0.5D0
-      R(0,3)=cos(30.D0*deg2rad)
-      R(1,3)=0.D0
-      R(2,3)=0.D0
-      R(0,4)=0.D0
-      R(1,4)=0.D0
-      R(2,4)=0.5D0
-      R(0,5)=0.D0
-      R(1,5)=cos(26.D0*deg2rad)
-      R(2,5)=0.5D0+sin(26.D0*deg2rad)
-      do i=1,5
-        do j=0,2
-          R(j,i)=vbl*R(j,i)
-        enddo
-      enddo
-c      i=move_res(R(0,1),0.D0*deg2rad,180.D0*deg2rad)
-      imov=2
-      i=move_res(0.D0*deg2rad,180.D0*deg2rad,imov)
-      print *,'RETURNED ',i
-      print *,(R(i,3)/vbl,i=0,2)
-
-      return
-      end
-
-c-------------------------------------------------------------