anchor.diff:            x     = x + xp*Srefl
             y     = y + yp*Srefl
   
  -          Dxp = (((-three*Ang_rms)*xp_rel)/(two*xpcrit) + Ang_avr) + ((three*Ang_rms)*(L_chan/c_rcurv))/(two*xpcrit)
  +          Dxp = (((-three*Ang_rms)*xp_rel)/(two*xpcrit) + Ang_avr) + ((three*Ang_rms)*(tdefl))/(two*xpcrit)
             xp  = xp + Dxp
             x   = x + (half*xp)*(s_length-Srefl)
             y   = y + (half*yp)*(s_length-Srefl)
  
positive.code:       ! of the particle -not along the longitudinal coordinate...
        if(ldech < l_chan) then
          iProc = proc_DC
          Dxp   = Ldech/c_rcurv ! Change angle from channeling [mrad]
          Sdech = Ldech*cos_mb(c_miscut + half*Dxp)
          x     = x  + Ldech*(sin_mb(half*Dxp+c_miscut)) ! Trajectory at channeling exit
          xp    = xp + Dxp + (two*(rndm4()-half))*xpcrit
  
negative.code:   z  = z - zp*s_rot
    s  = zero
  
    write(*,*) 'x after reference change', x, xp
  
    ! Check that particle hit the crystal
    if(x >= zero .and. x < c_xmax) then
  
  