reorganized +polytens (sorry manas)

matlab/+polytens/+cheb/coefs2vals.m

@@ -0,0 +1,4 @@
+function bmat = coefs2vals(norder,uvs)
+    bmat = polytens.cheb.pols(norder,uvs);
+    bmat = bmat.';
+end

matlab/+polytens/cheb_ders.m → matlab/+polytens/+cheb/ders.m

@@ -1,4 +1,4 @@
-function [p,dersu,dersv] = cheb_pols(norder,uv)
+function [p,dersu,dersv] = ders(norder,uv)
    [~,m] = size(uv);
    chebpols_u = zeros(m,norder+1);
    chebpols_v = zeros(m,norder+1);

matlab/+polytens/cheb_nodes.m → matlab/+polytens/+cheb/nodes.m

@@ -1,4 +1,4 @@
-function xvals = cheb_nodes(norder,kind)
+function xvals = nodes(norder,kind)
     kind0 = 1;
     if(nargin == 2)
         kind0 = kind;

matlab/+polytens/cheb_pols.m → matlab/+polytens/+cheb/pols.m

@@ -1,4 +1,4 @@
-function p = cheb_pols(norder,uv)
+function p = pols(norder,uv)
    [~,m] = size(uv);
    chebpols_u = zeros(m,norder+1);
    chebpols_v = zeros(m,norder+1);

matlab/+polytens/+cheb/vals2coefs.m

@@ -0,0 +1,4 @@
+function amat = vals2coefs(norder,uvs)
+    bmat = polytens.cheb.coefs2vals(norder,uvs);
+    amat = inv(bmat);
+end

matlab/+polytens/cheb_weights.m → matlab/+polytens/+cheb/weights.m

@@ -1,4 +1,4 @@
-function wts = cheb_weights(norder)
+function wts = weights(norder)
     [~,w] = polytens.cheb.pts(norder+1);
     wts = w'*w;
     wts = wts(:).'; 

matlab/+polytens/+lege/coefs2vals.m

@@ -0,0 +1,4 @@
+function bmat = coefs2vals(norder,uvs)
+    bmat = polytens.lege.pols(norder,uvs);
+    bmat = bmat.';
+end

matlab/+polytens/lege_ders.m → matlab/+polytens/+lege/ders.m

@@ -1,4 +1,4 @@
-function [p, dersu, dersv] = lege_ders(norder,uv)
+function [p, dersu, dersv] = ders(norder,uv)
    [~,m] = size(uv);
 
    legpols_u = zeros(m,norder+1);

matlab/+polytens/lege_nodes.m → matlab/+polytens/+lege/nodes.m

@@ -1,4 +1,4 @@
-function xvals = lege_nodes(norder)
+function xvals = nodes(norder)
     [x] = polytens.lege.pts(norder+1);
     [yy,xx] = meshgrid(x);
     xx = xx(:);

matlab/+polytens/lege_pols.m → matlab/+polytens/+lege/pols.m

@@ -1,4 +1,4 @@
-function p = lege_pols(norder,uv)
+function p = pols(norder,uv)
    [~,m] = size(uv);
    legpols_u = zeros(m,norder+1);
    legpols_v = zeros(m,norder+1);

matlab/+polytens/+lege/vals2coefs.m

@@ -0,0 +1,4 @@
+function amat = vals2coefs(norder,uvs)
+    bmat = polytens.lege.coefs2vals(norder,uvs);
+    amat = inv(bmat);
+end

matlab/+polytens/lege_weights.m → matlab/+polytens/+lege/weights.m

@@ -1,4 +1,4 @@
-function wts = lege_weights(norder)
+function wts = weights(norder)
     [~,w] = polytens.lege.pts(norder+1);
     wts = w'*w;
     wts = wts(:).'; 

matlab/@surfer/interpolate_data.m

@@ -0,0 +1,66 @@
+function dens_int = interpolate_data(obj, dens, ipatch_ids, uvs_targ)
+% INTERPOLATE_DATA interpolates a density on surfer object 
+%
+% dens_int = interpolate_data(obj,dens,ipatch_id,uvs_targ) where dens_int is a 
+% vector with values of dens interpolated to the given ipatch_ids and uvs_targ
+%
+% Input: obj    - a surfer object
+%        dens   - a density defined on obj
+%        ipatch_id - patch_ids of targs
+%        uvs_targ - uv_coordinates
+%
+% Output: values of the interpolated density
+
+    if size(dens,1) == obj.npts
+       fl = false;
+    elseif size(dens,2) == obj.npts
+       dens = dens.';
+       fl = true;
+    else
+        error("INTERPOLATE_DATA:error in size");
+    end
+
+    upids = unique(ipatch_ids);
+    npatch = length(upids);
+    dens_int = zeros([length(ipatch_ids), size(dens,2)]);
+
+for ii = 1:npatch
+
+    patch_id = upids(ii);
+    norder = obj.norders(patch_id);
+    iptype = obj.iptype(patch_id);
+
+    vals = dens(obj.patch_id == patch_id,:);
+    uvs = obj.uvs_targ(:,obj.patch_id == patch_id);
+    uvs_int = uvs_targ(:,ipatch_ids == patch_id);
+
+    if iptype == 1
+     %   Rokhlin-Vioreanu nodes
+
+     v2c = koorn.vals2coefs(norder, uvs);
+     c2v = koorn.coefs2vals(norder, uvs_int);
+     dens_int(ipatch_ids == patch_id,:) = c2v*(v2c*vals);
+
+    elseif iptype == 11
+     %   Gauss-Legendre nodes
+
+    v2c = polytens.lege.vals2coefs(norder,uvs);
+    c2v = polytens.lege.coefs2vals(norder,uvs_int);
+    dens_int(ipatch_ids == patch_id,:) = c2v*(v2c*vals);
+
+    elseif iptype == 12
+     %   Chebyshev nodes
+
+    v2c = polytens.cheb.vals2coefs(norder,uvs);
+    c2v = polytens.cheb.coefs2vals(norder,uvs_int);
+    dens_int(ipatch_ids == patch_id,:) = c2v*(v2c*vals);
+
+    end
+
+end
+
+if fl
+    dens_int = dens_int .';
+end
+
+end

matlab/@surfer/oversample.m

@@ -64,15 +64,15 @@
 
 
         elseif(iptype == 11)
-            rnodes = polytens.lege_nodes(norder);
-            rnodes_over = polytens.lege_nodes(nover);
-            vmats{i} = polytens.lege_coefs2vals(norder,rnodes_over);
-            umats{i} = polytens.lege_vals2coefs(norder,rnodes);
+            rnodes = polytens.lege.nodes(norder);
+            rnodes_over = polytens.lege.nodes(nover);
+            vmats{i} = polytens.lege.coefs2vals(norder,rnodes_over);
+            umats{i} = polytens.lege.vals2coefs(norder,rnodes);
         elseif(iptype == 12)
-            rnodes = polytens.cheb_nodes(norder);
-            rnodes_over = polytens.cheb_nodes(nover);
-            vmats{i} = polytens.cheb_coefs2vals(norder,rnodes_over);
-            umats{i} = polytens.cheb_vals2coefs(norder,rnodes);
+            rnodes = polytens.cheb.nodes(norder);
+            rnodes_over = polytens.cheb.nodes(nover);
+            vmats{i} = polytens.cheb.coefs2vals(norder,rnodes_over);
+            umats{i} = polytens.cheb.vals2coefs(norder,rnodes);
         end
         if(norder>nover)
             [lia,iindb] = ismembertol(rnodes_over',rnodes',1e-7,'ByRows',true);

matlab/@surfer/plot.m

@@ -95,15 +95,15 @@
         rnodes = zeros(2,nuse);
         rnodes(1,:) = x(:);
         rnodes(2,:) = y(:);
-        pols{i} = polytens.lege_pols(norder, rnodes);
+        pols{i} = polytens.lege.pols(norder, rnodes);
         nstot{i} = nuse;
         T{i} = tri_quad;
     elseif ip0 == 12
         nuse = length(x(:));
         rnodes = zeros(2,nuse);
         rnodes(1,:) = x(:);
         rnodes(2,:) = y(:);
-        pols{i} = polytens.cheb_pols(norder, rnodes);
+        pols{i} = polytens.cheb.pols(norder, rnodes);
         nstot{i} = nuse;
         T{i} = tri_quad;
     end

matlab/@surfer/surfer.m

@@ -59,6 +59,8 @@
 %                              representation
 %   rotate(obj, eul)         - rotate surface based on euler angles
 %   scale(obj, sf)           - scale object
+%   interpolate_data(obj,...)- interpolate density to uv_targs and
+%                              patch_ids
 %   translate(obj, r)        - translate object
 %   merge([array of objs])   - merge an array of surface objects
 %   area(obj)                - compute the surface area of object
@@ -183,21 +185,21 @@
                     epts = [0,0,1;0,1,0];
                     p_vert_mats{i} = koorn.pols(no_ip_uni(i,1), epts);
                 elseif(ip0 == 11)
-                    rnodes{i} = polytens.lege_nodes(no_ip_uni(i,1));
-                    rwts{i} = polytens.lege_weights(no_ip_uni(i,1));
-                    umats{i} = polytens.lege_vals2coefs(no_ip_uni(i,1),rnodes{i});
-                    [~, dumats{i}, dvmats{i}] = polytens.lege_ders(no_ip_uni(i,1), rnodes{i});
+                    rnodes{i} = polytens.lege.nodes(no_ip_uni(i,1));
+                    rwts{i} = polytens.lege.weights(no_ip_uni(i,1));
+                    umats{i} = polytens.lege.vals2coefs(no_ip_uni(i,1),rnodes{i});
+                    [~, dumats{i}, dvmats{i}] = polytens.lege.ders(no_ip_uni(i,1), rnodes{i});
                     npols{i} = size(rnodes{i},2);
                     epts = [-1,1,1,-1;-1,-1,1,1];
-                    p_vert_mats{i} = polytens.lege_pols(no_ip_uni(i,1), epts);
+                    p_vert_mats{i} = polytens.lege.pols(no_ip_uni(i,1), epts);
                 elseif(ip0 == 12)
-                    rnodes{i} = polytens.cheb_nodes(no_ip_uni(i,1));
-                    rwts{i} = polytens.cheb_weights(no_ip_uni(i,1));
-                    umats{i} = polytens.cheb_vals2coefs(no_ip_uni(i,1),rnodes{i});
-                    [~, dumats{i}, dvmats{i}] = polytens.cheb_ders(no_ip_uni(i,1), rnodes{i});
+                    rnodes{i} = polytens.cheb.nodes(no_ip_uni(i,1));
+                    rwts{i} = polytens.cheb.weights(no_ip_uni(i,1));
+                    umats{i} = polytens.cheb.vals2coefs(no_ip_uni(i,1),rnodes{i});
+                    [~, dumats{i}, dvmats{i}] = polytens.cheb.ders(no_ip_uni(i,1), rnodes{i});
                     npols{i} = size(rnodes{i},2);
                     epts = [-1,1,1,-1;-1,-1,1,1];
-                    p_vert_mats{i} = polytens.cheb_pols(no_ip_uni(i,1), epts);
+                    p_vert_mats{i} = polytens.cheb.pols(no_ip_uni(i,1), epts);
                 else
                     fprintf('Invalid type of patch, exiting\n');
                     return
@@ -324,6 +326,7 @@
          a = area(obj);
          [srcvals,srccoefs,norders,ixyzs,iptype,wts] = extract_arrays(obj);
          [objout,varargout] = oversample(obj,novers);
+         dens_int = interpolate_data(obj, dens, ipatch_ids, uvs_targ);
          [objout,varargout] = affine_transf(obj,mat,shift);
          [varargout] = scatter(obj,varargin);
          [spmat] = conv_rsc_to_spmat(obj,row_ptr,col_ind,wnear);

matlab/@surfer/vals2coefs.m

@@ -52,11 +52,11 @@
             rnodes = koorn.rv_nodes(norder);
             umats{i} = koorn.vals2coefs(norder,rnodes);       
         elseif(iptype == 11)
-            rnodes = polytens.lege_nodes(norder);
-            umats{i} = polytens.lege_vals2coefs(norder,rnodes);
+            rnodes = polytens.lege.nodes(norder);
+            umats{i} = polytens.lege.vals2coefs(norder,rnodes);
         elseif(iptype == 12)
-            rnodes = polytens.cheb_nodes(norder);
-            umats{i} = polytens.cheb_vals2coefs(norder,rnodes);
+            rnodes = polytens.cheb.nodes(norder);
+            umats{i} = polytens.cheb.vals2coefs(norder,rnodes);
         end
         npp  = size(umats{i},2);
         inds = (intmp == i);

matlab/src/vals_to_coefs_surface_fun.m

@@ -34,11 +34,11 @@
             rnodes = koorn.rv_nodes(norder);
             umats{i} = koorn.vals2coefs(norder,rnodes);  
         elseif(iptype == 11)
-            rnodes = polytens.lege_nodes(norder);
-            umats{i} = polytens.lege_vals2coefs(norder,rnodes);
+            rnodes = polytens.lege.nodes(norder);
+            umats{i} = polytens.lege.vals2coefs(norder,rnodes);
         elseif(iptype == 12)
-            rnodes = polytens.cheb_nodes(norder);
-            umats{i} = polytens.cheb_vals2coefs(norder,rnodes);
+            rnodes = polytens.cheb.nodes(norder);
+            umats{i} = polytens.cheb.vals2coefs(norder,rnodes);
         end
     end
 

matlab/test/test_interpolate_data.m

@@ -0,0 +1,57 @@
+test_interpolate_data0();
+
+function test_interpolate_data0()
+
+iseed = 319;
+rng(iseed);
+
+S = geometries.disk([],[],[3 3 3],8,1);
+dens = [eval_gauss(S.r(1,:),S.r(2,:)) eval_gauss(S.r(2,:),S.r(1,:)).^2];
+
+ipatchids = randsample(S.npatches,4);
+uvs_targ = [0.1 0.8; 0.3 0.5; 0.7 0.1; 0.24 0.31].';
+
+xs = S.interpolate_data(S.r(1,:),ipatchids,uvs_targ);
+ys = S.interpolate_data(S.r(2,:),ipatchids,uvs_targ);
+vals = S.interpolate_data(dens,ipatchids,uvs_targ);
+
+vals2 = [eval_gauss(xs,ys) eval_gauss(xs,ys).^2];
+
+assert(vecnorm( vals(:) - vals2(:) ) < 1e-9)
+
+S = geometries.disk([],[],[3 3 3],8,11);
+
+iseed = 312;
+rng(iseed);
+dens = [eval_gauss(S.r(1,:),S.r(2,:)) eval_gauss(S.r(2,:),S.r(1,:)).^2].';
+
+ipatchids = randsample(S.npatches,4);
+
+xs = S.interpolate_data(S.r(1,:),ipatchids,uvs_targ);
+ys = S.interpolate_data(S.r(2,:),ipatchids,uvs_targ);
+vals = S.interpolate_data(dens,ipatchids,uvs_targ);
+
+vals2 = [eval_gauss(xs,ys) eval_gauss(xs,ys).^2].';
+
+assert(vecnorm( vals(:) - vals2(:) ) < 1e-9)
+
+S = geometries.disk([],[],[3 3 3],8,12);
+
+iseed = 312;
+rng(iseed);
+dens = [eval_gauss(S.r(1,:),S.r(2,:)) eval_gauss(S.r(2,:),S.r(1,:)).^2];
+
+xs = S.interpolate_data(S.r(1,:),ipatchids,uvs_targ);
+ys = S.interpolate_data(S.r(2,:),ipatchids,uvs_targ);
+vals = S.interpolate_data(dens,ipatchids,uvs_targ);
+
+vals2 = [eval_gauss(xs,ys) eval_gauss(xs,ys).^2];
+
+assert(vecnorm( vals(:) - vals2(:) ) < 1e-9)
+
+
+end
+
+function val = eval_gauss(x,y)
+    val = exp(-x(:).^2 - y(:).^2);
+end