anchor.diff:              Dxp = (Length-Lrefl)/c_rcurv
               x   = x + sin_mb(half*Dxp+xp)*Rlength ! Trajectory at channeling exit
               y   = y + red_S*yp
  -            xp  = Length/c_rcurv + (half*ran_gauss(one))*xpcrit ! [mrad]
  +            xp  = Length/c_rcurv + (half*ran_gauss(zero))*xpcrit ! [mrad]
   
               call cry_calcIonLoss(is,pc,rlength,dest)
               pc = pc - (half*dest)*Rlength  ! "added" energy loss once captured
  
positive.code:       teta = sqrt(t)/pc
      end if
  
      tx = (teta*ran_gauss(one))*c1e3
      tz = (teta*ran_gauss(one))*c1e3
  
      ! Change p angle
      xp = xp + tx
  
negative.code:           call cry_calcIonLoss(is,pc,length,dest)
            pc = pc - dest*length ! energy loss to ionization [GeV]
          else
            Dxp = L_chan/c_rcurv + (half*ran_gauss(zero))*xpcrit ! Change angle[rad]
            xp  = Dxp
            x   = x + L_chan*(sin_mb(half*Dxp+c_miscut)) ! Trajectory at channeling exit
            y   = y + s_length * yp
  
            call cry_calcIonLoss(is,pc,length,dest)
  