anchor.diff:      sigmaX  = scatter_proList(idPro)%fParams(12)
       sigmaY  = scatter_proList(idPro)%fParams(13)
       normFac = scatter_proList(idPro)%fParams(14)
  -    ! sRot    = scatter_proList(idPro)%fParams(15)
  -    ! cRot    = scatter_proList(idPro)%fParams(16)
  -    ! sOrbXP  = scatter_proList(idPro)%fParams(18)
  -    ! sOrbYP  = scatter_proList(idPro)%fParams(19)
  -
  -    ! xRot    = (xPos - orbX + sPos*sOrbXP)*cRot - (yPos - orbY + sPos*sOrbXP)*sRot
  -    ! yRot    = (xPos - orbX + sPos*sOrbYP)*sRot + (yPos - orbY + sPos*sOrbYP)*cRot
  +
       xProj   = xPos - orbX
       yProj   = yPos - orbY
  -    ! xRot    = xProj*cRot - yProj*sRot
  -    ! yRot    = xProj*sRot + yProj*cRot
  -    rhoOut  = (nBeam*normFac) * exp_mb(-half*(xProj/sigmaX)**2 -half*(yProj/sigmaY)**2)
  +    rhoOut  = (nBeam*normFac) * exp_mb(    &
  +      - half*(xProj/sigmaX)**2             &
  +      - half*(yProj/sigmaY)**2             &
  +      - half*(sPos/scatter_beam2Length)**2 &
  +    )
   
     case default
   
  
positive.code:     ! orbY   = scatter_proList(idPro)%fParams(9)
  
        ! Compute the sigmas from beam 2 emittance
        sigmaX    = sqrt((scatter_beam2EmitX * betaX) / (gamma0 / beta0)) * c1m3
        sigmaY    = sqrt((scatter_beam2EmitY * betaY) / (gamma0 / beta0)) * c1m3
        sigmaXeff = sqrt(sigmaX**2 + (scatter_beam2Length * tan_mb(orbXP * c1m3))**2)
        sigmaYeff = sqrt(sigmaY**2 + (scatter_beam2Length * tan_mb(orbYP * c1m3))**2)
  
        ! Compute approximate effective sigmas from crossing angle, assuming sigma_z >> sigma_x,y
        ! xFac  = sqrt(one + ((scatter_beam2Length / sigmaX) * tan_mb(orbXP * c1m3))**2)
  
negative.code:       ! sigmaYeff = abs(sigmaX*cRot - sigmaY*sRot)
  
        ! Compute approximate kinematic factor and the normalisation of the PDF
        kFac    = two*cos_mb(two*(sigmaXeff/scatter_beam2Length))*cos_mb(two*(sigmaYeff/scatter_beam2Length))
        normFac = kFac/(twopi*sqrt(sigmaXeff**2 * sigmaYeff**2))
  
        scatter_proList(idPro)%fParams(10) = sigmaX
        scatter_proList(idPro)%fParams(11) = sigmaY
  