Looks like your data is virtually free of noise. That, combined with a very high sampling frequency would mean that at the jumps the data is exactly at the threshold between two quantized values. Set up nodes at the middle points of the vertical jumps and construct splines that connect the nodes. The easiest is to just draw straight lines between successive nodes, which gives a piece-wise constant differential. I wonder if that is good enough, or if you already tried that. If you need the velocity real-time, this approach is problematic because occasionally you might have to wait for a new node for some time.
You can further low-pass filter the interpolated data. If you use a filter with an impulse response that is nowhere negative, such as a Gaussian function, then there will be no overshoot.
With linear interpolation, everywhere between successive nodes, the speed will be simply the position difference between the nodes divided by the time difference between the nodes. You can run the smoothing filter on that piece-wise constant speed data and the result will be the same as if you'd run it on the linear position ramps and then differentiated (associative property of convolution, as both differentiation and filtering are convolution).