1 | % initialize data with ideal solution
|
---|
2 | x0=.05
|
---|
3 | y0=.05
|
---|
4 | [x,y]=meshgrid(linspace(-1,1,10), linspace(-1,1,10))
|
---|
5 | dv=(x(1,2)-x(1,1)).^2;
|
---|
6 | mass_target=pi*(1-exp(-1))
|
---|
7 | rho_predict=arrayfun(@(x,y) exp(-((x-x0)^2+(y-y0)^2)).*((x-x0).^2+(y-y0).^2<1), x, y)/mass_target;
|
---|
8 | mask=((x-x0).^2+(y-y0).^2 < 1);
|
---|
9 |
|
---|
10 |
|
---|
11 | % Plot initial guess
|
---|
12 | subplot(1,4,1)
|
---|
13 | imagesc(x(1,:),y(:,1),rho_predict)
|
---|
14 | colorbar()
|
---|
15 |
|
---|
16 | %Predict correct mass, and center
|
---|
17 | mass_predict=sum(sum(rho_predict))*dv
|
---|
18 | com_x=sum(sum(rho_predict.*x))*dv/mass_predict
|
---|
19 | com_y=sum(sum(rho_predict.*y))*dv/mass_predict
|
---|
20 | title(['initial guess - com = (',num2str(com_x),',',num2str(com_y),')', ' mass=', num2str(mass_predict)])
|
---|
21 |
|
---|
22 |
|
---|
23 | % Calculate solution in least squares sense
|
---|
24 | A=[x(mask)';y(mask)';ones(1,numel(x(mask)))]*dv;
|
---|
25 | B=[x0; y0; mass_target];
|
---|
26 | sol=pinv(A)*B;
|
---|
27 | rho_correct=zeros(size(rho_predict));
|
---|
28 | rho_correct(mask)=sol;
|
---|
29 |
|
---|
30 | subplot(1,4,2)
|
---|
31 | imagesc(x(1,:), y(:,1), rho_correct)
|
---|
32 | colorbar()
|
---|
33 | title('Solution with smallest norm')
|
---|
34 |
|
---|
35 | sol=lsqlin(eye(numel(x(mask))), rho_predict(mask), [],[],A,B)
|
---|
36 | subplot(1,4,3)
|
---|
37 | rho_correct=zeros(size(rho_predict));
|
---|
38 | rho_correct(mask)=sol;
|
---|
39 | imagesc(x(1,:), y(:,1),rho_correct)
|
---|
40 | title('Solution closest to prediction')
|
---|
41 | colorbar()
|
---|
42 |
|
---|
43 |
|
---|
44 | sol=rho_predict(mask)-A'*inv(A*A')*(A*rho_predict(mask) - B);
|
---|
45 | subplot(1,4,4)
|
---|
46 | rho_correct=zeros(size(rho_predict));
|
---|
47 | rho_correct(mask)=sol;
|
---|
48 | imagesc(x(1,:), y(:,1),rho_correct)
|
---|
49 | title('Solution closest to prediction')
|
---|
50 | colorbar()
|
---|