Honestly, I don't see the problem, at least if the discontinuity is as large as you show it. That will make it clear where the breaks live. So then your model will be just a a base curve, plus a couple of Heaviside functions with estimated coefficients. I'll make up some data...
X = sort(rand(30,1))*1.5;
H = @(X,origin) +(X >= origin);
Y = sin(X) + H(X,Bloc(1))*Bval(1) + H(X,Bloc(2))*Bval(2) + randn(size(X))/40;
So a curve with two breaks in it. We can never know where the break occurs, except that it MUST occur between two points. So just search for the two largest jumps in value between any pair of points. Then we can assume the breaks happened at the mid point between the located points. Since I've sorted the X values, we can just use diff.
[~,Jind] = maxk(abs(diff(Y)),2);
Blocest = sort((X(Jind) + X(Jind+1))/2)
So not too bad. Again, we cannot know exactly where that break really occurs, but we can assume it falls midway between the indicated points.
But now we can estimate a model. For example, a cubic polynomial will probably be adequate for a curve this simple.
A = [X(:).^(polydeg:-1:0),H(X(:),Blocest(1)),H(X(:),Blocest(2))];
So a cubic polynomial, plus the two discontinuity terms. See it did a pretty good job of estimating the jumps.
Ypred = polyval(coef(1:end-2),X) + coef(end-1)*H(X,Blocest(1)) + coef(end)*H(X,Blocest(2));
It seems ok to me.
The curve, with breaks removed, is just:
Ycorrected = Y - coef(end-1)*H(X,Blocest(1)) - coef(end)*H(X,Blocest(2));
If your data is far more noisy than you show, or the discontinuities are far smaller, then you may have a problem. But then you may just be better off fitting a curve with perhaps only one jump discontinuity. If you can't see a jump, then is it a jump at all? You might use a tolerance to decide how significant the jumps must be for it to recognized as such.
Finally, if the width of the jogged segment is known in advance, then you can even get trickier about how you find the locations of the breaks.