Function to undertake steps 2 to 5 of the GDR algorithm (Clause 8)
function [at, bt, da, db, t, xs, z, f, g, h, F2, g0, h0, gt, ht, Gt2, r] = ... algm_gdr2_steps_2_to_5(x, ux, y, uy, uxy, at, bt, da, db, ... t, xs, z, f, g, h, F2, g0, h0, gt, ht, Gt2, r, ind)
Step 2.
m = length(x); t{ind} = ones(m, 1)./(uy.^2 -2*bt{ind}*uxy + bt{ind}^2*ux.^2); xs{ind} = (x.*(uy.^2 - bt{ind}*uxy) - (uxy - bt{ind}*ux.^2).*(y - at{ind})).*t{ind}; z{ind} = y - at{ind} - bt{ind}*x;
Step 3.
f{ind} = sqrt(t{ind}); g{ind} = f{ind}.*xs{ind}; h{ind} = f{ind}.*z{ind};
Step 4.
Substep (i).
F2{ind} = sum(f{ind}.*f{ind});
Substep (ii).
g0{ind} = sum(f{ind}.*g{ind})/F2{ind}; h0{ind} = sum(f{ind}.*h{ind})/F2{ind};
Substep (iii).
gt{ind} = g{ind} - g0{ind}*f{ind}; ht{ind} = h{ind} - h0{ind}*f{ind};
Substep (iv).
Gt2{ind} = sum(gt{ind}.*gt{ind});
Substep (v).
db{ind} = sum(gt{ind}.*ht{ind})/Gt2{ind}; da{ind} = h0{ind} - db{ind}*g0{ind};
Step 5.
at{ind+1} = at{ind} + da{ind}; bt{ind+1} = bt{ind} + db{ind}; r{ind} = ht{ind} - db{ind}*gt{ind};
End of algm_gdr2_steps_2_to_5.m