anchor.diff:    x_i = x_i - Np*dP    ! Rescale the incoming x at the left crystalline plane
     x_i = x_i - (dP/two) ! Rescale the incoming x in the middle of crystalline planes
   
  -  pv   = pc*c1e9*pc*c1e9/sqrt(pc*c1e9*pc*c1e9 + 93828.0e6_fPrec*93828.0e6_fPrec) ! Calculate pv=P/E
  -  Ueff = eUm(is)*(two*x_i/dp)*(two*x_i/dp)+pv*x_i/r                              ! Calculate effective potential
  -  Et   = pv*xp*xp/two+Ueff                                                       ! Calculate transverse energy
  -  Ec   = eUm(is)*(one-rc/r)*(one-rc/r)                                           ! Calculate critical energy in bent crystals
  +  pv   = pc**2/sqrt(pc**2 + (pmap*c1m3)**2)*c1e9    ! Calculate pv=P/E
  +  Ueff = eUm(is)*(two*x_i/dp)*(two*x_i/dp)+pv*x_i/r ! Calculate effective potential
  +  Et   = pv*xp*xp/two+Ueff                          ! Calculate transverse energy
  +  Ec   = eUm(is)*(one-rc/r)*(one-rc/r)              ! Calculate critical energy in bent crystals
   
     ! To avoid negative Et
     xminU = -dp*dp*pc*c1e9/(eight*eUm(is)*r)
  
positive.code: 
    real(kind=fPrec), parameter :: re  = 2.818d-15   ! electron radius [m]
    real(kind=fPrec), parameter :: me  = 0.510998910 ! electron mass [MeV/c^2]
    real(kind=fPrec), parameter :: mp  = 938.272013  ! proton mass [MeV/c^2]
  
    real(kind=fPrec), parameter :: aTF = 0.194d-10   ! Screening function [m]
    real(kind=fPrec), parameter :: dP  = 1.92d-10    ! Distance between planes (110) [m]
  
negative.code:   use coll_common, only : coll_debug
    use coll_crystal
    use coll_materials, only : nmat, rho, anuc, hcut, bnref, csref, csect
    use physical_constants
  
    implicit none
  
  