% measure some box geometries. code from talkfig_orderNstrip.m
% Barnett 9/25/12
clear;
N = 130;  % base # pts
s = segment.smoothnonsym(N, 0.3, 0.3, 3);
d = domain(s, 1); k = 20.0; d.k = k;  % choose Laplace or Helmholtz
o = []; d.addlayerpot(s, 'D', o); l = d.bas{1};

delta = 2.5;  % how many "h" away center is from boundary
al = delta/N; % imag strip width wrt [0,1]-scaled seg
Bdelta = 2.5; B = ceil(N/(2*Bdelta)); % how many boxes
t0 = (1:B)/B+1i*al; z0 = s.Z(t0);  % all preimages and centers

if 1, gam = nan(B,1); figure; s.plot; hold on;  % plot geom ........
for b=1:B
  box = l.QBXbox(real(t0(b)), [1/(2*B) 1/(1*B)]); % note not always 2x1 boxes
  R = max(abs(box.v - z0(b))); % max radius of box about z0
  gam(b) = s.QBXgammaconst(t0(b), R);
  plot([box.v; box.v(1)], 'b-'); plot(z0(b), 'k*');
end
hold on; plot(z0(find(gam==max(gam))), 'ro'); % worst box
figure; plot(1:B, gam, '+-'); title('gamma consts for z0''s');
median(gam), max(gam)
end % ..........

dx = 0.01; bb = d.boundingbox;  % choose eval grid...
[xx yy] = meshgrid(bb(1):dx:bb(2),bb(3):dx:bb(4));
zz = xx(:)+1i*yy(:); tic; di = d.inside(zz); ii = find(di);
fprintf('inside %d pts, %d verts = %.3g s\n',numel(zz),numel(s.approxv),toc);
