anchor.diff:        s_int = (-b_eq - sqrt(delta))/(two*a_eq)
         s_imp = s_int
   
  +      ! MISCUT first step: P coordinates
  +      s_P_tmp = (c_rcurv-c_xmax)*sin_mb(-c_miscut)
  +      x_P_tmp = c_xmax + (c_rcurv-c_xmax)*cos_mb(-c_miscut)
  +      write(*,*) "P side tmp", s_P_tmp, x_P_tmp
  +
         if(s_int < c_rcurv*c_sBend) then
           ! Transform to a new reference system: shift and rotate
           x_int  = xp*s_int + x
  
positive.code:   real(kind=fPrec) Pvr                  ! Prob for VR->AM transition
  
    real(kind=fPrec) const_dech,xpin,ypin,tchan,L_chan,mep
  
    real(kind=fPrec), parameter :: c_v1 =  1.7_fPrec ! Fitting coefficient
    real(kind=fPrec), parameter :: c_v2 = -1.5_fPrec ! Fitting coefficient
  
negative.code:   real(kind=fPrec), intent(inout) :: yp
    real(kind=fPrec), intent(inout) :: pc
    real(kind=fPrec), intent(in)    :: length
    real(kind=fPrec), intent(in)    :: s_P
    real(kind=fPrec), intent(in)    :: x_P
  
    integer nam,zn                        ! Switch on/off the nuclear interaction (NAM) and the MCS (ZN)
    real(kind=fPrec) ymax,ymin            ! Crystal geometrical parameters
  