source: trunk/src/tps_eval.m @ 363

Last change on this file since 363 was 356, checked in by sommeria, 13 years ago

civ updated with new functions for opening files, consistently with uvmat
Bugs to be expected (use previous version then)

File size: 1.4 KB
Line 
1%'tps_eval': calculate the thin plate spline (tps) interpolation at a set of points
2% see tps_ceff for more information
3%------------------------------------------------------------------------
4% function EM = tps_eval(dsites,ctrs)
5%------------------------------------------------------------------------
6% OUPUT:
7% EM:  Mx(N+s) matrix representing the contributions at the M sites
8%   from unit sources located at each of the N centers, + (s+1) columns
9%   representing the contribution of the linear gradient part.
10%
11%INPUT:
12%dsites:  Nxs matrix representing the postions of the N 'observation' sites, with s the space dimension
13%ctrs: Mxs matrix  representing the postions of the M centers, sources of the tps,
14%
15% related functions:
16% tps_coeff, tps_eval_dxy
17
18function EM = tps_eval(dsites,ctrs)
19[M,s] = size(dsites); [N,s] = size(ctrs);
20EM = zeros(M,N);
21
22% calculate distance matrix: accumulate sum of squares of coordinate differences
23% The ndgrid command produces two MxN matrices:
24%   Dsite, consisting of N identical columns (each containing
25%       the d-th coordinate of the M data sites)
26%   Ctrs, consisting of M identical rows (each containing
27%       the d-th coordinate of the N centers)
28for d=1:s
29 [Dsites,Ctrs] = ndgrid(dsites(:,d),ctrs(:,d));
30 EM = EM + (Dsites-Ctrs).^2;%EM=square of distance matrices
31end
32
33% calculate tps
34np=find(EM~=0);
35EM(np) = EM(np).*log(EM(np))/2;%= tps formula r^2 log(r) (EM=r^2)
36
37% add linear gradient part:
38EM = [EM ones(M,1) dsites];
Note: See TracBrowser for help on using the repository browser.