I'm going to make an assumption that the functions have the same domain and are meant to represent the same thing.
My suggested approach would be: for each datapoint in the original function, find its corresponding datapoint on the smoothed spline. Find the vector between each point (if this is just a 1D equation [f(x)=y] then it will just be the differences in y) and take its magnitude. Then, you could sum the magnitudes of the vectors in order to quanitfy how close or far the step function is from the smooth spline. If the step function is continuous, it may be necessary to discretize it so that it is of the form f below.
E.g.
% Our related dataset
f = [(0:pi/20:pi/2)', rand(11,1)];
% Our smoothed function
g = [(0:pi/20:pi/2)', sin((0:pi/20:pi/2))'];
% Our summated difference function based on the dataset and smoothed
% function. We desire a zero deviation from smoothness; the higher the
% number, the less smooth the step function is relative to the spline
devFromSmoothness = sum(sqrt(sum((f-g).^2,2)),1);
Cheers,
Jakeb
