Create a function handle that is a sum of function handles

Goodmorning everyone,

I have the following points

POINTS=[0.0194358127458200	0.0125146090698240	0.0120324235527160	0.00761174198214000	0.00660775699279100	0.00381485098617300	0.00195617023840500	0.000824508358196000	0.00121681256884900	0.000390389885216000	0.000444315057222000	0;
                         0	0.000452782917043000	0.000395689955442000	0.00121172812659200	0.000829869449880000	0.00195546440425700	0.00381482017361200	0.00661064156533300	0.00762496776011900	0.0120411026638820	0.0125098828784870	0.0189137987850580];

I need to create a function handle as follows.

I know that the center of the circumference that approximates better those points is function of R, and is found in the point [R,R]. The approximation error is then:

error=@(R) abs(norm([R-POINTS(1,i),R-POINTS(2,i)])-R));

This for every point.

I want to create the function handle that is the sum of all the errors like:

1) totdist=@(R) error1+error2+error3+...

I can create a single function function manually by summing every error like:

2) totdist=@(R) abs(norm([R-POINTS(1,1),R-POINTS(2,1)])-R))+abs(norm([R-POINTS(1,2),R-POINTS(2,2)])-R))+abs(norm([R-POINTS(1,3),R-POINTS(2,3)])-R))+...

However when the number of points become large, I would like to create the summing function 1) authomatically without manually summing 100 times as in 2).

The function must be used with "fminsearch" afterwards.

Is there a way to achieve this?

Thanks in advance for the answer.

[Bjorn Gustavsson]

@Riccardo Mari: at least try to work through the suggestions to see what they do. Use simple test-cases.

Answer 1

[Bjorn Gustavsson]

Your specification is a bit ambiguous to me, but if I interpret it such that you want to fit a function such that:

(I've used the additional condition that and , but this one can be easily relaxed). You can make a general solution something like this:

f_x = @(R,theta,f_theta) R + f_theta(R,theta).*cos(theta);
f_y = @(R,theta,f_theta) R + f_theta(R,theta).*sin(theta);
err_fcn = @(R,X,Y,f_theta) sum( (X-f_x(R,atan2(Y-R,X-R),f_theta)).^2 + (Y-f_y(R,atan2(Y-R,X-R),f_theta)).^2);
f_theta = @(R,theta) R*ones(size(theta));
Rest = fminsearch(@(R) err_fcn(R,POINTS(1,:)',POINTS(2,:)',f_theta),20e-3);

This is a rather cumbersome way to make a simple circle-fitting-function, you will find multiple circle-fitting-functions on the file exchange: circle-fitting FEx contributions. If your circumference is a bit more complicated, even more complicated than an ellipse even, then you simply have to modify the f_theta function-handle accordingly. If the centre of the circumference is not exactly the same as the radius for example then you'll have to modify the above snippet a bit:

f_x = @(R,theta,f_theta) R(1) + f_theta(R(3),theta).*cos(theta);
f_y = @(R,theta,f_theta) R(2) + f_theta(R(3),theta).*sin(theta);
err_fcn = @(R,X,Y,f_theta) sum( (X-f_x(R,atan2(Y-R(2),X-R(1)),f_theta)).^2 + (Y-f_y(R,atan2(Y-R(2),X-R(1)),f_theta)).^2);
f_theta = @(R,theta) R*ones(size(theta));
Rest3 = fminsearch(@(R) err_fcn(R,POINTS(1,:)',POINTS(2,:)',f_theta),20e-3*[1 1 1]);

You can also allow for additional parameters to fit to, for example major and minor radius and angle of semi-major axis of an ellipse. But, this starts to become rather cluttered code for anonymous function-handles at this stage and somewhere here I tend to transition to proper matlab-functions in files, just to keep the code cleaner.

HTH

[Riccardo Mari]

I try to explain better.

I have to measure the error in the distance as

error=@(R) abs(norm([R-POINTS(1,i),R-POINTS(2,i)])-R));
% This is the difference between the distance of the points to the centre
% of the circle (in the point [R,R]) and the radius R.

So I have to apply this to every point in the list. R is the radius of the circle. I want to find the radius that minimizes the total error. So I have

error1=@(R) abs(norm([R-POINTS(1,1),R-POINTS(2,1)])-R));
error2=@(R) abs(norm([R-POINTS(1,2),R-POINTS(2,2)])-R));
error3=@(R) abs(norm([R-POINTS(1,3),R-POINTS(2,3)])-R));
%And so on until the last one. In this case:
error12=@(R) abs(norm([R-POINTS(1,12),R-POINTS(2,12)])-R));

I have to use "fminsearch" with the sum of those errors. For twelve points, like in this case, I have to manually write:

totdist=@(R) (abs(norm([R-POINTS(1,1),R-POINTS(2,1)])-R)+abs(norm([R-POINTS(1,2),R-POINTS(2,2)])-R)+abs(norm([R-POINTS(1,3),R-POINTS(2,3)])-R)+abs(norm([R-POINTS(1,4),R-POINTS(2,4)])-R)+abs(norm([R-POINTS(1,5),R-POINTS(2,5)])-R)+abs(norm([R-POINTS(1,6),R-POINTS(2,6)])-R)+abs(norm([R-POINTS(1,7),R-POINTS(2,7)])-R)+abs(norm([R-POINTS(1,8),R-POINTS(2,8)])-R)+abs(norm([R-POINTS(1,9),R-POINTS(2,9)])-R)+abs(norm([R-POINTS(1,10),R-POINTS(2,10)])-R)+abs(norm([R-POINTS(1,11),R-POINTS(2,11)])-R)+abs(norm([R-POINTS(1,12),R-POINTS(2,12)])-R));

You see that is quite long to manually write, especially when the points number becomes large.

My question is if there is a way to generate all the error_n functions and then sum them to get the function handle totdist without writing manually every function, like using a for loop to generate each error function and the sum function to sum them together.

Hope it's clearer now.

[Bjorn Gustavsson]

That's what my first example does. For general circle-fitting functions follow the link to the good file-exchange contributions that solve that problem for the more general case.

HTH

Answer 2

[Matt J]

error=@(R) abs( vecnorm(POINTS-R) - R );

[Matt J]

I have to measure the error in the distance as error=@(R) abs(norm([R-POINTS(1,i),R-POINTS(2,i)])-R)) .... So I have to apply this to every point in the list.

No, you do not have to implement a separate error function for every Point(:,i). There are vectorized functions which calculate all errors at once, which is what my answer has given you.

I want to find the radius that minimizes the total error.

If you know the center [x0,y0] of the circle already, the minimizing R is just the average of the empricial distances. There is no need for fminsearch.

R = mean(vecnorm(POINTS-[x0;y0]))

Answer 3

[Riccardo Mari]

Thanks everybody for the answers.

I'm sorry if I couldn't be clearer because I can't share too much abount the code and the project.

I used the solution proprosed here. To summarize, I did this:

N=length(POINTS);
f=cell(N,1);
for i=1:N
    f{i}=@(R) abs(norm([R-POINTS(1,i),R-POINTS(2,i)])-R);
end
totdist=@(R) summer(f,R);

with:

function y=summer(f,R)
y=0;
for i=1:numel(f)
    y=y+f{i}(R);
end
end

[Matt J]

It's a bad solution. Very slow and inefficient. I urge you to reconsider the other proposals.

POINTS=rand(2,100000);  R=0;

tic;
    N=length(POINTS);
    f=cell(N,1);
    for i=1:N
        f{i}=@(R) abs(norm([R-POINTS(1,i),R-POINTS(2,i)])-R);
    end
    totdist=@(R) summer(f,R);
    
    
    version1=totdist(R); %calculation
toc;

Elapsed time is 0.279596 seconds.

tic
    errorFcn=@(R) norm( vecnorm(POINTS-R) - R ,1);
    
    version2=errorFcn(R);%calculation
toc

Elapsed time is 0.006358 seconds.

version1,version2   %same result

version1 = 7.6622e+04

version2 = 7.6622e+04