diff --git a/chunkie/+chnk/+flam/kernbyindexr.m b/chunkie/+chnk/+flam/kernbyindexr.m
index 8f6b1eb..a03fa2a 100644
--- a/chunkie/+chnk/+flam/kernbyindexr.m
+++ b/chunkie/+chnk/+flam/kernbyindexr.m
@@ -1,4 +1,4 @@
-function mat = kernbyindexr(i,j,targs,chnkr,kern,opdims,spmat)
+function mat = kernbyindexr(i,j,targobj,chnkr,kern,opdims,spmat)
 %% evaluate system matrix by entry index utility function for 
 % general kernels, with replacement for specific entries and the
 % ability to add a low-rank modification, rectangular version
@@ -17,6 +17,10 @@
 %
 % i - array of row indices to compute
 % j - array of col indices to compute
+% targobj - object describing the target points, can be specified as
+%       * array of points
+%       * chunker object
+%       * chunkgraph object
 % chnkr - chunker object describing boundary
 % whts - smooth integration weights on chnkr
 % kern - kernel function of the form kern(s,t,stau,ttau) where s and t 
@@ -42,15 +46,26 @@
 [juni,~,ijuni] = unique(jpts);
 
 % matrix-valued entries of kernel for unique points
-
-ri = targs(:,iuni); rj = chnkr.r(:,juni);
+rj = chnkr.r(:,juni);
 dj = chnkr.d(:,juni); d2j = chnkr.d2(:,juni);
 nj = chnkr.n(:,juni);
-%di = bsxfun(@rdivide,di,sqrt(sum(di.^2,1)));
-%dj = bsxfun(@rdivide,dj,sqrt(sum(dj.^2,1)));
 srcinfo = []; srcinfo.r = rj; srcinfo.d = dj; srcinfo.d2 = d2j;
 srcinfo.n = nj;
-targinfo = []; targinfo.r = ri;
+% Assign appropriate object to targinfo
+targinfo = [];
+if isa(targobj, "chunker")
+    targinfo.r = targobj.r(:,iuni);
+    targinfo.d = targobj.d(:,iuni);
+    targinfo.d2 = targobj.d2(:,iuni);
+    targinfo.n = targobj.n(:,iuni);
+elseif isa(targobj, "chunkgraph")
+    targinfo.r = targobj.r(:,iuni);
+    targinfo.d = targobj.d(:,iuni);
+    targinfo.d2 = targobj.d2(:,iuni);
+    targinfo.n = targobj.n(:,iuni);
+else
+    targinfo.r = targobj(:,iuni);
+end
 
 matuni = kern(srcinfo,targinfo);
 
diff --git a/chunkie/chunkerinterior.m b/chunkie/chunkerinterior.m
index 9563985..49ece7e 100644
--- a/chunkie/chunkerinterior.m
+++ b/chunkie/chunkerinterior.m
@@ -1,16 +1,18 @@
-function [in] = chunkerinterior(chnkr,pts,opts)
+function [in] = chunkerinterior(chnkobj,ptsobj,opts)
 %CHUNKERINTERIOR returns an array indicating whether each point specified 
 % % by pts is inside the domain. Assumes the domain is closed.
 %
-% Syntax: in = chunkerinterior(chnkr,pts,opts)
-%         in = chunkerinterior(chnkr,{x,y},opts) % meshgrid version
+% Syntax: in = chunkerinterior(chnkobj,pts,opts)
+%         in = chunkerinterior(chnkobj,{x,y},opts) % meshgrid version
 %
 % Input:
-%   chnkr - chunker object describing geometry
-%   pts - (chnkr.dim,:) array of points to test
-% 
-%   {x,y} - length 2 cell array. the points checked then have the
-%       coordinates of a mesh grid [xx,yy] = meshgrid(x,y)
+%   chnkobj - chunker object or chunkgraph object describing geometry
+%   ptsobj - object describing the target points, can be specified as
+%       * (chnkr.dim,:) array of points to test
+%       * {x,y} - length 2 cell array. the points checked then have the
+%           coordinates of a mesh grid [xx,yy] = meshgrid(x,y)
+%       * chunker object, in which case it uses chunker.r(:,:) 
+%       * chunkgraph object, in which case it uses chunkgraph.r(:,:) 
 %
 % Optional input:
 %   opts - options structure with entries:
@@ -34,23 +36,39 @@
 
 grid = false;
 
-assert(chnkr.dim == 2,'interior only well-defined for 2D');
 
 if nargin < 3
     opts = [];
 end
 
-if isa(pts,"cell")
-    assert(length(pts)==2,'second input should be either 2xnpts array or length 2 cell array');
-    x = pts{1};
-    y = pts{2};
-    grid = true;
+% Assign appropriate object to chnkr
+if class(chnkobj) == "chunker"
+   chnkr = chnkobj;
+elseif class(chnkobj) == "chunkgraph"
+   chnkr = merge(chnkobj.echnks);
+else
+    msg = "Unsupported object in chunkerinteriort";
+    error(msg)
 end
 
-if grid
+assert(chnkr.dim == 2,'interior only well-defined for 2D');
+
+% Assign appropriate object to pts
+if isa(ptsobj, "cell")
+    assert(length(ptsobj)==2,'second input should be either 2xnpts array or length 2 cell array');
+    x = ptsobj{1};
+    y = ptsobj{2};
+    grid = true;
     [xx,yy] = meshgrid(x,y);
     pts = [xx(:).'; yy(:).'];
-end    
+elseif isa(ptsobj, "chunker")
+    pts = ptsobj.r(:,:);
+elseif isa(ptsobj, "chunkgraph")
+    pts = ptsobj.r(:,:);
+else
+    pts = ptsobj;
+end
+
 
 usefmm = true;
 if isfield(opts,'fmm')
diff --git a/chunkie/chunkerkerneval.m b/chunkie/chunkerkerneval.m
index 24d068c..573bf92 100644
--- a/chunkie/chunkerkerneval.m
+++ b/chunkie/chunkerkerneval.m
@@ -1,17 +1,20 @@
-function fints = chunkerkerneval(chnkr,kern,dens,targs,opts)
+function fints = chunkerkerneval(chnkobj,kern,dens,targobj,opts)
 %CHUNKERKERNEVAL compute the convolution of the integral kernel with
 % the density defined on the chunk geometry. 
 %
 % Syntax: fints = chunkerkerneval(chnkr,kern,dens,targs,opts)
 %
 % Input:
-%   chnkr - chunker object description of curve
+%   chnkobj - chunker object or chunkgraph object description of curve
 %   kern - kernel class object or kernel function taking inputs 
 %                      kern(srcinfo,targinfo) 
 %   dens - density on boundary, should have size opdims(2) x k x nch
 %          where k = chnkr.k, nch = chnkr.nch, where opdims is the 
 %           size of kern for a single src,targ pair
-%   targs - targ(1:2,i) gives the coords of the ith target
+%   targobj - object describing the target points, can be specified as
+%       * array of points
+%       * chunker object
+%       * chunkgraph object
 %
 % Optional input:
 %   opts - structure for setting various parameters
@@ -54,6 +57,19 @@
     kerneval = kern;
 end
 
+% Assign appropriate object to chnkr
+if class(chnkobj) == "chunker"
+   chnkr = chnkobj;
+elseif class(chnkobj) == "chunkgraph"
+   chnkr = merge(chnkobj.echnks);
+else
+    msg = "Unsupported object in chunkerinteriort";
+    error(msg)
+end
+
+
+
+
 srcinfo = []; targinfo = [];
 srcinfo.r = chnkr.r(:,1); srcinfo.d = chnkr.d(:,1); 
 srcinfo.n = chnkr.n(:,1); srcinfo.d2 = chnkr.d2(:,1);
@@ -86,31 +102,49 @@
 if isfield(opts,'fac'); opts_use.fac = opts.fac; end
 if isfield(opts,'eps'); opts_use.eps = opts.eps; end
 
-[dim,~] = size(targs);
+% Assign appropriate object to targinfo
+targinfo = [];
+if isa(targobj, "chunker")
+    targinfo.r = targobj.r(:,:);
+    targinfo.d = targobj.d(:,:);
+    targinfo.d2 = targobj.d2(:,:);
+    targinfo.n = targobj.n(:,:);
+elseif isa(targobj, "chunkgraph")
+    targinfo.r = targobj.r(:,:);
+    targinfo.d = targobj.d(:,:);
+    targinfo.d2 = targobj.d2(:,:);
+    targinfo.n = targobj.n(:,:);
+else
+    targinfo.r = targobj;
+end
+
+
+
+[dim,~] = size(targinfo.r);
 
 if (dim ~= 2); warning('only dimension two tested'); end
 
 if opts_use.forcesmooth
-    fints = chunkerkerneval_smooth(chnkr,kern,opdims,dens,targs, ...
+    fints = chunkerkerneval_smooth(chnkr,kern,opdims,dens,targinfo, ...
         [],opts_use);
     return
 end
 
 if opts_use.forceadap
     fints = chunkerkerneval_adap(chnkr,kern,opdims,dens, ...
-        targs,[],opts_use);
+        targinfo,[],opts_use);
     return
 end
 
 
 if opts_use.forcepquad
     optsflag = []; optsflag.fac = opts_use.fac;
-    flag = flagnear(chnkr,targs,optsflag);
-    fints = chunkerkerneval_smooth(chnkr,kern,opdims,dens,targs, ...
+    flag = flagnear(chnkr,targinfo.r,optsflag);
+    fints = chunkerkerneval_smooth(chnkr,kern,opdims,dens,targinfo, ...
         flag,opts_use);
 
     fints = fints + chunkerkerneval_ho(chnkr,kern,opdims,dens, ...
-        targs,flag,opts_use);
+        targinfo,flag,opts_use);
 
     return
 end
@@ -120,27 +154,26 @@
 %optsflag = []; optsflag.fac = opts_use.fac;
 rho = 1.8;
 optsflag = [];  optsflag.rho = rho;
-flag = flagnear_rectangle(chnkr,targs,optsflag);
+flag = flagnear_rectangle(chnkr,targinfo.r,optsflag);
 
 npoly = chnkr.k*2;
 nlegnew = chnk.ellipse_oversample(rho,npoly,opts_use.eps);
 nlegnew = max(nlegnew,chnkr.k);
 
 [chnkr2,dens2] = upsample(chnkr,nlegnew,dens);
-fints = chunkerkerneval_smooth(chnkr2,kern,opdims,dens2,targs, ...
+fints = chunkerkerneval_smooth(chnkr2,kern,opdims,dens2,targinfo, ...
     flag,opts_use);
 
 fints = fints + chunkerkerneval_adap(chnkr,kern,opdims,dens, ...
-        targs,flag,opts_use);
-
+        targinfo,flag,opts_use);
 
 end
 
 function fints = chunkerkerneval_ho(chnkr,kern,opdims,dens, ...
-    targs,flag,opts)
+    targinfo,flag,opts)
 
 % target
-[~,nt] = size(targs);
+[~,nt] = size(targinfo.r);
 fints = zeros(opdims(1)*nt,1);
 k = size(chnkr.r,2);
 [t,w] = lege.exps(2*k);
@@ -155,8 +188,22 @@
 % interpolation matrix
 intp = lege.matrin(k,t);          % interpolation from k to 2*k
 intp_ab = lege.matrin(k,[-1;1]);  % interpolation from k to end points
+
+targs = targinfo.r;
+
 targd = zeros(chnkr.dim,nt); targd2 = zeros(chnkr.dim,nt);
 targn = zeros(chnkr.dim,nt);
+if isfield(targinfo, 'd')
+    targd = targinfo.d;
+end
+
+if isfield(targinfo, 'd2')
+    targd2 = targinfo.d2;
+end
+
+if isfield(targinfo, 'n')
+    targn = targinfo.n;
+end
 for j=1:size(chnkr.r,3)
     [ji] = find(flag(:,j));
     if(~isempty(ji))
@@ -174,7 +221,7 @@
 
 
 function fints = chunkerkerneval_smooth(chnkr,kern,opdims,dens, ...
-    targs,flag,opts)
+    targinfo,flag,opts)
 
 if isa(kern,'kernel')
     kerneval = kern.eval;
@@ -203,12 +250,10 @@
 dens = reshape(dens,opdims(2),k,nch);
 
 [~,w] = lege.exps(k);
-[~,nt] = size(targs);
+[~,nt] = size(targinfo.r);
 
 fints = zeros(opdims(1)*nt,1);
 
-targinfo = []; targinfo.r = targs;
-
 % assume smooth weights are good enough
 
 % Sequence of checks, first see ifflam is set as it supercedes
@@ -272,24 +317,6 @@
     end
 else
 
-%     % hack to ignore interactions: create sparse array with very small
-%     % number instead of zero in ignored entries. kernbyindexr overwrites
-%     %   with that small number
-%     [i,j] = find(flag);
-%     
-%     % targets correspond to multiple outputs (if matrix valued kernel)
-%     inew = (opdims(1)*(i(:)-1)).' + (1:opdims(1)).'; inew = inew(:);
-%     jnew = (j(:)).' + 0*(1:opdims(1)).'; jnew = jnew(:);
-%     
-%     % source chunks have multiple points and can be multi-dimensional
-%     jnew = (opdims(2)*chnkr.k*(jnew(:)-1)).' + (1:(chnkr.k*opdims(2))).';
-%     jnew = jnew(:);
-%     inew = (inew(:)).' + 0*(1:(chnkr.k*opdims(2))).';
-%     inew = inew(:);
-%     
-%     mm = nt*opdims(1); nn = chnkr.npt*opdims(2);
-%     v = 1e-300*ones(length(inew),1); sp = sparse(inew,jnew,v,mm,nn);
-
     wts = chnkr.wts;
     wts = wts(:);
     
@@ -297,14 +324,29 @@
         xflam1 = chnkr.r(:,:);
         xflam1 = repmat(xflam1,opdims(2),1);
         xflam1 = reshape(xflam1,chnkr.dim,numel(xflam1)/chnkr.dim);
-        targsflam = repmat(targs(:,:),opdims(1),1);
-        targsflam = reshape(targsflam,chnkr.dim,numel(targsflam)/chnkr.dim);
-        matfun = @(i,j) chnk.flam.kernbyindexr(i,j,targs,chnkr,kerneval, ...
-            opdims);
+
+        targinfo_flam = [];
+        targinfo_flam.r = repelem(targinfo.r(:,:),1,opdims(1));
+        if isfield(targinfo, 'd')
+            targinfo_flam.d = repelem(targinfo.d(:,:),1,opdims(1));
+        end
+        
+        if isfield(targinfo, 'd2')
+            targinfo_flam.d2 = repelem(targinfo.d2(:,:),1,opdims(1));
+        end
+        
+        if isfield(targinfo, 'n')
+            targinfo_flam.n = repelem(targinfo.n(:,:),1,opdims(1));
+        end
+
+% TODO: Pull through data?
+
+        matfun = @(i,j) chnk.flam.kernbyindexr(i, j, targinfo_flam, ...,
+                           chnkr, kerneval, opdims);
     
 
         width = max(abs(max(chnkr)-min(chnkr)))/3;
-        tmax = max(targs(:,:),[],2); tmin = min(targs(:,:),[],2);
+        tmax = max(targinfo.r(:,:),[],2); tmin = min(targinfo.r(:,:),[],2);
         wmax = max(abs(tmax-tmin));
         width = max(width,wmax/3);  
         npxy = chnk.flam.nproxy_square(kerneval,width);
@@ -314,16 +356,16 @@
             ctr,chnkr,kerneval,opdims,pr,ptau,pw,pin);
 
         optsifmm=[]; optsifmm.Tmax=Inf;
-        F = ifmm(matfun,targsflam,xflam1,200,1e-14,pxyfun,optsifmm);
+        F = ifmm(matfun,targinfo_flam.r,xflam1,200,1e-14,pxyfun,optsifmm);
         fints = ifmm_mv(F,dens(:),matfun);
     else
         wts2 = repmat(wts(:).', opdims(2), 1);
         sigma = wts2(:).*dens(:);
-        fints = kern.fmm(1e-14, chnkr, targs(:,:), sigma);
+        fints = kern.fmm(1e-14, chnkr, targinfo.r(:,:), sigma);
     end
     % delete interactions in flag array (possibly unstable approach)
     
-    targinfo = [];
+    
     if ~isempty(flag)
         for i = 1:nch
             densvals = dens(:,:,i); densvals = densvals(:);
@@ -338,9 +380,24 @@
             delsmooth = find(flag(:,i)); 
             delsmoothrow = (opdims(1)*(delsmooth(:)-1)).' + (1:opdims(1)).';
             delsmoothrow = delsmoothrow(:);
-            targinfo.r = targs(:,delsmooth);
 
-            kernmat = kerneval(srcinfo,targinfo);
+            targinfo_use = [];
+            targinfo_use.r = targinfo.r(:,delsmooth);
+
+            if isfield(targinfo, 'd')
+                targinfo_use.d = targinfo.d(:,delsmooth);
+            end
+        
+            if isfield(targinfo, 'd2')
+                targinfo_use.d2 = targinfo.d2(:,delsmooth);
+            end
+        
+            if isfield(targinfo, 'n')
+                targinfo_use.n = targinfo.n(:,delsmooth);
+            end
+
+
+            kernmat = kerneval(srcinfo,targinfo_use);
             fints(delsmoothrow) = fints(delsmoothrow) - kernmat*densvals;
         end
     end    
@@ -349,7 +406,7 @@
 end
 
 function fints = chunkerkerneval_adap(chnkr,kern,opdims,dens, ...
-    targs,flag,opts)
+    targinfo,flag,opts)
 
 if isa(kern,'kernel')
     kerneval = kern.eval;
@@ -370,7 +427,23 @@
 assert(numel(dens) == opdims(2)*k*nch,'dens not of appropriate size')
 dens = reshape(dens,opdims(2),k,nch);
 
+% Extract target info
+targs = targinfo.r;
 [~,nt] = size(targs);
+targd = zeros(chnkr.dim,nt); targd2 = zeros(chnkr.dim,nt);
+targn = zeros(chnkr.dim,nt);
+if isfield(targinfo, 'd')
+    targd = targinfo.d;
+end
+
+if isfield(targinfo, 'd2')
+    targd2 = targinfo.d2;
+end
+
+if isfield(targinfo, 'n')
+    targn = targinfo.n;
+end
+
 
 fints = zeros(opdims(1)*nt,1);
 
@@ -384,8 +457,8 @@
 n = chnkr.n;
 d2 = chnkr.d2;
 h = chnkr.h;
-targd = zeros(chnkr.dim,nt); targd2 = zeros(chnkr.dim,nt);
-targn = zeros(chnkr.dim,nt);
+
+
 
 if isempty(flag) % do all to all adaptive
     for i = 1:nch
diff --git a/chunkie/chunkerkernevalmat.m b/chunkie/chunkerkernevalmat.m
index c0a96c9..31b9526 100644
--- a/chunkie/chunkerkernevalmat.m
+++ b/chunkie/chunkerkernevalmat.m
@@ -1,4 +1,4 @@
-function mat = chunkerkernevalmat(chnkr,kern,targs,opts)
+function mat = chunkerkernevalmat(chnkr,kern,targobj,opts)
 %CHUNKERKERNEVALMAT compute the matrix which maps density values on 
 % the chunk geometry to the value of the convolution of the given
 % integral kernel with the density at the specified target points
@@ -8,7 +8,10 @@
 % Input:
 %   chnkr - chunker object description of curve
 %   kern - integral kernel taking inputs kern(srcinfo,targinfo) 
-%   targs - targ(1:2,i) gives the coords of the ith target
+%   targobj - object describing the target points, can be specified as
+%       * array of points
+%       * chunker object
+%       * chunkgraph object
 %
 % Optional input:
 %   opts - structure for setting various parameters
@@ -68,7 +71,23 @@
 if isfield(opts,'fac'); fac = opts.fac; end
 if isfield(opts,'eps'); eps = opts.eps; end
 
-[dim,~] = size(targs);
+% Assign appropriate object to targinfo
+targinfo = [];
+if isa(targobj, "chunker")
+    targinfo.r = targobj.r(:,:);
+    targinfo.d = targobj.d(:,:);
+    targinfo.d2 = targobj.d2(:,:);
+    targinfo.n = targobj.n(:,:);
+elseif isa(targobj, "chunkgraph")
+    targinfo.r = targobj.r(:,:);
+    targinfo.d = targobj.d(:,:);
+    targinfo.d2 = targobj.d2(:,:);
+    targinfo.n = targobj.n(:,:);
+else
+    targinfo.r = targobj;
+end
+
+[dim,~] = size(targinfo.r);
 
 
 if (dim ~= 2); warning('only dimension two tested'); end
@@ -78,23 +97,23 @@
 optsadap.eps = eps;
 
 if forcesmooth
-    mat = chunkerkernevalmat_smooth(chnkr,kern,opdims,targs, ...
+    mat = chunkerkernevalmat_smooth(chnkr,kern,opdims,targinfo, ...
         [],optssmooth);
     return
 end
 
 if forceadap
     mat = chunkerkernevalmat_adap(chnkr,kern,opdims, ...
-        targs,[],optsadap);
+        targinfo,[],optsadap);
     return
 end
 
 if forcepqud
     optsflag = []; optsflag.fac = fac;
-    flag = flagnear(chnkr,targs,optsflag);
+    flag = flagnear(chnkr,targinfo.r,optsflag);
     spmat = chunkerkernevalmat_ho(chnkr,kern,opdims, ...
-        targs,flag,optsadap);
-    mat = chunkerkernevalmat_smooth(chnkr,kern,opdims,targs, ...
+        targinfo,flag,optsadap);
+    mat = chunkerkernevalmat_smooth(chnkr,kern,opdims,targinfo, ...
         flag,opts);
     mat = mat + spmat;
     return
@@ -103,16 +122,16 @@
 % smooth for sufficiently far, adaptive otherwise
 
 optsflag = []; optsflag.fac = fac;
-flag = flagnear(chnkr,targs,optsflag);
+flag = flagnear(chnkr,targinfo.r,optsflag);
 spmat = chunkerkernevalmat_adap(chnkr,kern,opdims, ...
-        targs,flag,optsadap);
+        targinfo,flag,optsadap);
 
 if nonsmoothonly
     mat = spmat;
     return;
 end
 
-mat = chunkerkernevalmat_smooth(chnkr,kern,opdims,targs, ...
+mat = chunkerkernevalmat_smooth(chnkr,kern,opdims,targinfo, ...
     flag,opts);
 
 mat = mat + spmat;
@@ -123,7 +142,7 @@
 
 
 function mat = chunkerkernevalmat_smooth(chnkr,kern,opdims, ...
-    targs,flag,opts)
+    targinfo,flag,opts)
 
 if nargin < 6
     flag = [];
@@ -135,7 +154,6 @@
 k = chnkr.k;
 nch = chnkr.nch;
 
-targinfo = []; targinfo.r = targs;
 srcinfo = []; srcinfo.r = chnkr.r(:,:); srcinfo.n = chnkr.n(:,:);
 srcinfo.d = chnkr.d(:,:); srcinfo.d2 = chnkr.d2(:,:);
 
@@ -163,7 +181,7 @@
 end
 
 function mat = chunkerkernevalmat_adap(chnkr,kern,opdims, ...
-    targs,flag,opts)
+    targinfo,flag,opts)
 
 k = chnkr.k;
 nch = chnkr.nch;
@@ -175,12 +193,30 @@
     opts = [];
 end
 
+% Extract target info
+targs = targinfo.r;
 [~,nt] = size(targs);
+targd = zeros(chnkr.dim,nt); targd2 = zeros(chnkr.dim,nt);
+targn = zeros(chnkr.dim,nt);
+if isfield(targinfo, 'd')
+    targd = targinfo.d;
+end
+
+if isfield(targinfo, 'd2')
+    targd2 = targinfo.d2;
+end
+
+if isfield(targinfo, 'n')
+    targn = targinfo.n;
+end
+
 
 % using adaptive quadrature
 
 
 if isempty(flag)
+
+
     mat = zeros(opdims(1)*nt,opdims(2)*chnkr.npt);
 
     [t,w] = lege.exps(2*k+1);
@@ -191,7 +227,7 @@
     n = chnkr.n;
     d2 = chnkr.d2;
     h = chnkr.h;
-    targd = zeros(chnkr.dim,nt); targd2 = zeros(chnkr.dim,nt);    
+    
     for i = 1:nch
         jmat = 1 + (i-1)*k*opdims(2);
         jmatend = i*k*opdims(2);
@@ -228,8 +264,6 @@
     n = chnkr.n;
     d2 = chnkr.d2;
     h = chnkr.h;
-    targd = zeros(chnkr.dim,nt); targd2 = zeros(chnkr.dim,nt);
-    targn = zeros(chnkr.dim,nt);
     for i = 1:nch
         jmat = 1 + (i-1)*k*opdims(2);
         jmatend = i*k*opdims(2);
@@ -260,7 +294,7 @@
 end
 
 function mat = chunkerkernevalmat_ho(chnkr,kern,opdims, ...
-    targs,flag,opts)
+    targinfo,flag,opts)
 
 k = chnkr.k;
 nch = chnkr.nch;
@@ -272,7 +306,24 @@
     opts = [];
 end
 
+% Extract target info
+targs = targinfo.r;
 [~,nt] = size(targs);
+targd = zeros(chnkr.dim,nt); targd2 = zeros(chnkr.dim,nt);
+targn = zeros(chnkr.dim,nt);
+if isfield(targinfo, 'd')
+    targd = targinfo.d;
+end
+
+if isfield(targinfo, 'd2')
+    targd2 = targinfo.d2;
+end
+
+if isfield(targinfo, 'n')
+    targn = targinfo.n;
+end
+
+
 
 % using Helsing-Ojala quadrature
 if isempty(flag) % figure out what is this flag for in adaptive routine
@@ -295,8 +346,7 @@
     % interpolation matrix 
     intp = lege.matrin(k,t);          % interpolation from k to 2*k
     intp_ab = lege.matrin(k,[-1;1]);  % interpolation from k to end points
-    targd = zeros(chnkr.dim,nt); targd2 = zeros(chnkr.dim,nt);
-    targn = zeros(chnkr.dim,nt);
+    
     for i = 1:nch
         jmat = 1 + (i-1)*k*opdims(2);
         jmatend = i*k*opdims(2);
diff --git a/chunkie/trapperkerneval.m b/chunkie/trapperkerneval.m
index d425b26..aaf7939 100644
--- a/chunkie/trapperkerneval.m
+++ b/chunkie/trapperkerneval.m
@@ -1,10 +1,21 @@
-function fints = trapperkerneval(trap,kern,dens,targs,opts)
+function fints = trapperkerneval(trap,kern,dens,targobj,opts)
 %TRAPPERINTKERN compute the convolution of the integral kernel with
 
 wts = trap.wts;
 srcinfo = []; srcinfo.r = trap.r; srcinfo.d = trap.d; srcinfo.d2 = trap.d2;
 srcinfo.n = trap.n;
-targinfo = []; targinfo.r = targs;
+
+% Assign appropriate object to targinfo
+targinfo = [];
+if isa(targobj, "trapper")
+    targinfo.r = targobj.r(:,:);
+    targinfo.d = targobj.d(:,:);
+    targinfo.d2 = targobj.d2(:,:);
+    targinfo.n = targobj.n(:,:);
+else
+    targinfo.r = targobj;
+end
+
 mat = kern(srcinfo,targinfo);
 fints = mat*diag(wts)*dens;
 
diff --git a/devtools/test/chunkerinteriorTest.m b/devtools/test/chunkerinteriorTest.m
index 65a6a91..2538be6 100644
--- a/devtools/test/chunkerinteriorTest.m
+++ b/devtools/test/chunkerinteriorTest.m
@@ -67,5 +67,20 @@
 opts = [];
 opts.fmm = true;
 opts.flam = false;
-start = tic; in3 = chunkerinterior(chnkr,{x1,x1},opts); t4 = toc(start);
-fprintf('%5.2e s : time for chunkerinterior (meshgrid, with fmm)\n',t4);
\ No newline at end of file
+start = tic; in4 = chunkerinterior(chnkr,{x1,x1},opts); t4 = toc(start);
+fprintf('%5.2e s : time for chunkerinterior (meshgrid, with fmm)\n',t4);
+
+% Test targets specified as chunker
+narms = 3;
+amp = 0.1;
+chnkr2 = chunkerfunc(@(t) 0.3*starfish(t,narms,amp),cparams,pref); 
+
+opts = [];
+opts.fmm = true;
+opts.flam = false;
+start = tic; in5 = chunkerinterior(chnkr,chnkr2,opts); t5 = toc(start);
+fprintf('%5.2e s : time for chunkerinterior (chunker, with fmm)\n',t5);
+assert(all(in5(:) == 1));
+
+
+
diff --git a/devtools/test/chunkgrphconstructTest.m b/devtools/test/chunkgrphconstructTest.m
index d829e77..6a6e50f 100644
--- a/devtools/test/chunkgrphconstructTest.m
+++ b/devtools/test/chunkgrphconstructTest.m
@@ -110,7 +110,7 @@
 
 prefs      = [];
 % prefs.chsmall = 1d-4;
-[cgrph] = chunkgraphinit(verts,edge2verts,fchnks,prefs);
+[cgrph] = chunkgraph(verts,edge2verts,fchnks,prefs);
 
 vstruc = procverts(cgrph);
 rgns = findregions(cgrph);
diff --git a/devtools/test/chunkgrphrcipTest.m b/devtools/test/chunkgrphrcipTest.m
index 154b123..391316a 100644
--- a/devtools/test/chunkgrphrcipTest.m
+++ b/devtools/test/chunkgrphrcipTest.m
@@ -62,7 +62,7 @@
 
 prefs      = [];
 % prefs.chsmall = 1d-4;
-[cgrph] = chunkgraphinit(verts,edge2verts,fchnks,prefs);
+[cgrph] = chunkgraph(verts,edge2verts,fchnks,prefs);
 
 vstruc = procverts(cgrph);
 rgns = findregions(cgrph);
@@ -110,7 +110,7 @@
 
 prefs      = [];
 % prefs.chsmall = 1d-4;
-[cgrph] = chunkgraphinit(verts,edge2verts,fchnks,prefs);
+[cgrph] = chunkgraph(verts,edge2verts,fchnks,prefs);
 
 vstruc = procverts(cgrph);
 rgns = findregions(cgrph);
diff --git a/devtools/test/chunkgrphrcipTransmissionTest.m b/devtools/test/chunkgrphrcipTransmissionTest.m
index 5547c62..7f1c4cf 100644
--- a/devtools/test/chunkgrphrcipTransmissionTest.m
+++ b/devtools/test/chunkgrphrcipTransmissionTest.m
@@ -55,7 +55,7 @@
 
 prefs      = [];
 % prefs.chsmall = 1d-4;
-[cgrph] = chunkgraphinit(verts,edge2verts,fchnks,prefs);
+[cgrph] = chunkgraph(verts,edge2verts,fchnks,prefs);
 
 
 vstruc = procverts(cgrph);
diff --git a/devtools/test/chunkrgrphOpdimTest.m b/devtools/test/chunkrgrphOpdimTest.m
index 4fb2652..38911c5 100644
--- a/devtools/test/chunkrgrphOpdimTest.m
+++ b/devtools/test/chunkrgrphOpdimTest.m
@@ -19,7 +19,7 @@
 fchnks = {}; % this by default gives me straight lines
 
 prefs = struct('maxchunklen',0.5);
-[cgrph] = chunkgraphinit(verts, edge2verts, fchnks, prefs);
+[cgrph] = chunkgraph(verts, edge2verts, fchnks, prefs);
 
 % figure(1); clf; hold on; axis equal; axis off;
 % plot(cgrph);
diff --git a/devtools/test/flamutilitiesTest.m b/devtools/test/flamutilitiesTest.m
index 0698b35..cd8d7bb 100644
--- a/devtools/test/flamutilitiesTest.m
+++ b/devtools/test/flamutilitiesTest.m
@@ -31,7 +31,7 @@
 end
 cparams = [];
 cparams.nover = 2;
-[cgrph] = chunkgraphinit(verts,edge2verts,fchnks,cparams);
+[cgrph] = chunkgraph(verts,edge2verts,fchnks,cparams);
 
 vstruc = procverts(cgrph);
 rgns = findregions(cgrph);