ectoSeg.m

% [im,imD] = MicroscopeData.Reader('H:/Smadar/Raw/22-06-2016_TimeLapse3_DMSO_20hrs_Calcein_FM464/matlab_decon', 'timeRange',[1,1]);
[im,imD] = MicroscopeData.Reader('C:/Users/mwinter/Documents/Lab/Smadar/Raw/22-06-2016_TimeLapse3_DMSO_20hrs_Calcein_FM464/matlab_decon', 'timeRange',[1,1]);

% addpath('C:\Users\mwinter\Documents\MATLAB\external\AOSLevelsetSegmentationToolboxM');
pxScale = imD.PixelPhysicalSize;

% mul = 20;
bwMask = llsm_mask(im);

se = strel('square', 31);
bwShrink = false(size(bwMask));
for z=1:size(bwMask,3)
    bwShrink(:,:,z) = imerode(bwMask(:,:,z), se);
end

imFM = mat2gray(im(:,:,:,1));

alpha = 0.95;
thresh = alpha*graythresh(imFM);
bwFM = imFM >= thresh;
bwFM(~bwShrink) = false;

% Approximate ectoderm thickness
ectoDist = 15.0;
idxPts = find(bwFM);
coords_xy = Utils.SwapXY_RC(Utils.IndToCoord(size(bwFM), idxPts));
coords_um = coords_xy.*pxScale;

shp = alphaShape(coords_um(:,1),coords_um(:,2),coords_um(:,3), ectoDist);
[f,pts_xyz] = boundaryFacets(shp);

unscalePts = pts_xyz ./ pxScale;

poly = D3d.Polygon.MakeEmptyStruct();
nrms = D3d.Polygon.CalcNorms(unscalePts, f);
poly(1).index = 1;
poly(1).frame = 1;
poly(1).label = '1';
poly(1).color = [1,0,0];
poly(1).faces = f;
poly(1).verts = unscalePts - 0.5;
poly(1).norms = nrms;
poly(1).CenterOfMass = mean(unscalePts);

D3d.Viewer.AddPolygons(poly);

% coords = Utils.IndToCoord(size(bwMask), find(bwMask));
% 
% init_phi = -ones(size(imFM));
% ctr = round(mean(coords));
% padding = 20;
% padctr = {ctr(1)-padding:ctr(1)+padding, ctr(2)-padding:ctr(2)+padding, ctr(3)-padding:ctr(3)+padding};
% init_phi(padctr{:}) = 1;
% 
% init_phi = ac_reinit(init_phi);
% 
% smooth_weight = 3;
% image_weight = 1e-3; 
% delta_t = 2;
% n_iters = 100;
% show_result = 1;
% 
% phi = ac_ChanVese_model(imFM, init_phi, smooth_weight, image_weight, delta_t, n_iters, show_result); 


% bwFM = mask;
% updatePolys(bwFM);
% 
% for i=1:500
%     bwFM = activecontour(imBMask, bwFM, 2, 'edge', 'ContractionBias', -20.0*mul, 'SmoothFactor',mul);
%     
%     updatePolys(bwFM);
% end

function updatePolys(bwFM)
    rp = regionprops(bwFM, 'PixelList');
    polys = [];
    for i=1:length(rp)
        polys = [polys; D3d.Polygon.Make(rp(i).PixelList, i, ['Endo_' num2str(i)], 1, [1 0 0], 1)];
    end

    D3d.Viewer.DeleteAllPolygons();
    D3d.Viewer.AddPolygons(polys);
    
    D3d.Update();
end