%% the saved file "*_res_cell.mat" contains the BW_cell
%% and cell centroid. In this one, the coordinates are stored in (Y,X) WITH the 0 in top-left
%% meaning that BW_cell(cell_centroid_imgRef(1,i), cell_centroid_imgRef(2,i), i) = 1 [counterintuitive]
function segment_cell
%% help
% comments
%% Parameters for positions and colors of figures and some handles
% The position of the game figure will be [figX1,figX2,figY1,figY2]
figX1 = 0;
figX2 = 1300;
figY1 = 0;
figY2 = 800;

%% sets initial variables
FileName = '';
PathName = '';
acquisition_name ='';
cell_centroid_imgRef = [];
frame = 1;
Nf = 150;
Img = [];
BW_cell = [];
openI = 1;
thresh = 0.001;


%% Create and hide the GUI figure as it is being constructed.
segfigure = figure('Visible','on','Tag','segfigure','Position',[figX1,figY1,figX2,figY2]);
set ( gcf, 'Color', [0 0 0] )

% To capture mouse and keyboard
set(segfigure,'WindowKeyPressFcn',@keypress);
set(segfigure,'WindowScrollWheelFcn',@scroll);
set(segfigure,'WindowButtonUpFcn',@click);
%% Create buttons and others
hmessages = uicontrol('Style','text','String','Please load a movie',...
    'Position',[50,730,350,50],'FontSize',14, ...
    'ForegroundColor','white','BackgroundColor','black'); % Display instructions

hNextSegbutton = uicontrol('Style','pushbutton',...
    'String','Segment Golgi','Callback',@hNextSeg_Callback,...
    'Position',[940,730,150,50],'FontSize',14); % Quit seg

hquitbutton = uicontrol('Style','pushbutton',...
    'String','Save and quit','Callback',@hquitbutton_Callback,...
    'Position',[1115,730,150,50],'FontSize',14); 

hloadbutton = uicontrol('Style','pushbutton',...
    'String','Load','Callback',@hloadbutton_Callback,...
    'Position',[400,730,100,50],'FontSize',14); % load file

hpath = uicontrol('Style','text','String','Path',...
    'Position',[510,730,400,50],'FontSize',10);

hframe = uicontrol('Style','slider',...
    'Min',1,'Max',Nf, 'value',frame,...
    'Position',[20,10,1200,20],...
    'SliderStep',[1/150 0.2],'Callback',@hframe_Callback);

hthresh = uicontrol('Style','slider',...
    'Min',-1,'Max',1, 'value',thresh,...
    'Position',[50,670,100,30],...
    'SliderStep',[0.4/200 0.1],'Callback',@hthresh_Callback);

hsegfow = uicontrol('Style','pushbutton','String', 'Seg >>',...
    'Position',[180,670,70,30],'FontSize',12,...
    'Callback',@hsegfow_Callback);

hautoth = uicontrol('Style','pushbutton','String', 'Auto th',...
    'Position',[280,670,70,30],'FontSize',12,...
    'Callback',@hautoth_Callback);


%% Create axes and display image
canvas = 600;
ImgOp = zeros(512,512);
sI = size(ImgOp);
panelI_X1 = 40;
panelI_X2 = canvas;
panelI_Y1 = 40;
panelI_Y2 = canvas;
ha = axes('Units','Pixels','Position',[panelI_X1,panelI_Y1,panelI_X2,panelI_Y2]);
hbw = axes('Units','Pixels','Position',[panelI_X1+canvas+30,panelI_Y1,panelI_X2,panelI_Y2]);

% Change units to normalized so components resize automatically.
set([segfigure,ha,hmessages...
    hbw,hpath,hloadbutton,hquitbutton,hNextSegbutton,...
    hframe,hsegfow...
    hthresh],...
    'Units','normalized');
%% Final settings
% Assign the GUI a name to appear in the window title.
set(segfigure,'Name','Segment Cell v 1.0')
% Move the GUI to the center of the screen.
movegui(segfigure,'center')
% Make the GUI visible.
set(segfigure,'Visible','on')

%% Callbacks
    function keypress(~,~)
        % Not used...
    end
    function click(~,~)
        % Not used...
    end
    function hNextSeg_Callback(~,~) % ajouter le lancement de l'autre GUI, enregistrer les contours de chaque cellule
        isOK_closing = checkSegmentation();
        
        if( isOK_closing == 1)
            computeResultsMatFile();
            close(segfigure); 
            
            save(fullfile(PathName,[acquisition_name 'res_cell.mat']),'cell_centroid_imgRef','BW_cell');
            segment_golgi(PathName,acquisition_name,BW_cell,cell_centroid_imgRef); 
        end
    end

    function segOk = checkSegmentation()
       segOk = 1;

        for i = 1:Nf
            bwp = bwperim(BW_cell(:,:,i));
            hold on
            [r,~] = find(bwp==1,1);
            if isempty(r) && segOk == 1
                frame = i;
                showI(Img(:,:,frame),BW_cell(:,:,frame));
                set(hmessages,'String',['Frame ',num2str(frame)]);
                set(hframe,'Value',frame);
                msgbox(['Image ' num2str(i) ' has no segmentation'],'Error','error');
                segOk = 0;
            end
        end
        
    end

    function computeResultsMatFile()
        h = waitbar(0,'Saving informations...');
        for framenumber = 1:Nf
            se = strel('disk',20);
            bwi = imclose(BW_cell(:,:,framenumber),se);
            bwi = imfill(bwi,'holes');
            
            stats = regionprops(im2bw(bwi),'Centroid','BoundingBox');
            waitbar(framenumber/Nf,h);
            if ~isempty(stats)
                ss = length(stats);
                if ss>1 % more than one reagions are found, the "most centered" is chosen
                    d=[];
                    for j=1:ss
                        d(j)=(stats(j).Centroid(1)-size(BW_cell,1)/2).^2+(stats(j).Centroid(2)-size(BW_cell,2)/2).^2;
                    end
                    [~,ig] = min(d);
                else
                    ig=1;
                end
                cell_centroid_imgRef(:,framenumber) = [stats(ig).Centroid(2),stats(ig).Centroid(1)];
                bw1reg = zeros(size(bwi));
                bw1reg(max(floor(stats(ig).BoundingBox(2)),1):min(ceil(stats(ig).BoundingBox(2)+stats(ig).BoundingBox(4)),size(bwi,1)), max(floor(stats(ig).BoundingBox(1)),1):min(ceil(stats(ig).BoundingBox(1)+stats(ig).BoundingBox(3)),size(bwi,2))) ...
                    = bwi(max(floor(stats(ig).BoundingBox(2)),1):min(ceil(stats(ig).BoundingBox(2)+stats(ig).BoundingBox(4)),size(bwi,1)), max(floor(stats(ig).BoundingBox(1)),1):min(ceil(stats(ig).BoundingBox(1)+stats(ig).BoundingBox(3)),size(bwi,2)));
                BW_cell(:,:,framenumber) = bw1reg;
            end
            
        end
        
        delete(h);
    end

    function hquitbutton_Callback(~,~) % ajouter le lancement de l'autre GUI, enregistrer les contours de chaque cellule
        isOK_closing = checkSegmentation();
        
        if( isOK_closing == 1)
           computeResultsMatFile();
            
            save(fullfile(PathName,[acquisition_name 'res_cell.mat']),'cell_centroid_imgRef','BW_cell');
            close(segfigure); 
        end
    end

    function hloadbutton_Callback(~,~)
        frame = 1;
        set(hframe,'Value',frame);
        [FileName, PathName] = uigetfile({'*.tif','*.tiff'}, 'Load the TIFF of your cell data');
        PathName = [PathName '\'];
        
        set(hpath,'String',PathName);
        
        if( ~isempty(PathName))
            tmp = strsplit(FileName,'_');
            acquisition_name = [tmp{1} '_' tmp{2} '_']; % res without extension and wavelength number

            info = imfinfo([PathName,FileName]);
            Nf = length(info);

            h = waitbar(0,'Please wait...');
            Img = ReadTiffStack([PathName,FileName]);

            BW_cell = zeros(size(Img,1),size(Img,2),size(Img,3));
            for i=1:Nf
                waitbar(i/Nf,h);
                BW_cell(:,:,i) = openi(Img(:,:,i));
                BW_cell(:,:,i) = segment(BW_cell(:,:,i));
            end

            delete(h);
            showI(Img(:,:,1),BW_cell(:,:,1));
            set(hframe,'max',Nf,'SliderStep',[1/(Nf-1) 0.2])
            set(hmessages,'String',['Loaded. Frame ',num2str(frame)]);
        end
    end

    function hframe_Callback(~,~)
        frame = get(hframe,'Value');
        frame = round(frame);
        set(hframe,'Value',frame);
        
        showI(Img(:,:,frame),BW_cell(:,:,frame));
        set(hmessages,'String',['Frame ',num2str(frame)]);
    end

    function hautoth_Callback(~,~)
        bwi = openi(Img(:,:,frame));
        bwi = bwi./max(max(bwi));
        th = graythresh(bwi);
        
        for i_t = -th:0.01:1-th
            bwo = im2bw(bwi,th+i_t);       
            bwo = imclearborder(bwo);
            bwo = bwareaopen(bwo,5000);
            se = strel('disk',10);
            bwo = imclose(bwo,se);
            bwo = imfill(bwo,'holes');
            if ( sum(bwo(:)) ~=0  )
                thresh = i_t;
                BW_cell(:,:,frame) = bwo;
            end
        end
        
        set(hthresh,'Value',thresh);
        showI(Img(:,:,frame),BW_cell(:,:,frame));
        
        if( sum(sum(BW_cell(:,:,frame))) ==0 )
            msgbox('Impossible to find a segmentation (try without imclearborder)','Error','error');
        end
        
    end

    function scroll(src,evnt)
        deltaframe = evnt.VerticalScrollCount;
        frame = frame + deltaframe;
        if frame<1
            frame=1;
        end
        if frame>Nf
            frame = Nf;
        end
        set(hframe,'Value',frame);
        showI(Img(:,:,frame),BW_cell(:,:,frame));
        set(hmessages,'String',['Frame ',num2str(frame)]);
    end

    function hthresh_Callback(~,~)
        thresh = get(hthresh,'Value');
        BW_cell(:,:,frame) = openi(Img(:,:,frame));
        BW_cell(:,:,frame) = segment(BW_cell(:,:,frame));
        axes(hbw);
        imshow(BW_cell(:,:,frame),[]);
        showI(Img(:,:,frame),BW_cell(:,:,frame));
    end

    function hsegfow_Callback(~,~)
        h = waitbar(0,'Please wait...');
        for i = frame:Nf
            waitbar(i/Nf,h)
            BW_cell(:,:,i) = openi(Img(:,:,i));
            BW_cell(:,:,i) = segment(BW_cell(:,:,i));
        end
        delete(h);
    end

%% Functions
    function showI(image,bwi)
        axes(ha);
        imshow(image,[]);
        hold on; 
        [~,contour] = computeCentroidOneImage(frame);
        plot(contour(:,2),contour(:,1),'g','LineWidth',0.5);
        hold off
        axes(hbw);
        imshow(bwi,[]);
        
    end


    function [bw1reg,cellcontour] = computeCentroidOneImage(framenumber)
        se = strel('disk',20);
        bwi = imclose(BW_cell(:,:,framenumber),se);
        bwi = imfill(bwi,'holes');
        
        stats = regionprops(im2bw(bwi),'Centroid','BoundingBox');
       
        if ~isempty(stats)
            ss = length(stats);
            
             if ss>1 % more than one reagions are found, the "most centered" is chosen
                    d=[];
                    for j=1:ss
                        d(j)=(stats(j).Centroid(1)-size(BW_cell,1)/2).^2+(stats(j).Centroid(2)-size(BW_cell,2)/2).^2;
                    end
                    [~,ig] = min(d);
                else
                    ig=1;
                end
           
            bw1reg = zeros(size(bwi));
            % keep only the bouynding box containing the 
            bw1reg(max(floor(stats(ig).BoundingBox(2)),1):min(ceil(stats(ig).BoundingBox(2)+stats(ig).BoundingBox(4)),size(bwi,1)), max(floor(stats(ig).BoundingBox(1)),1):min(ceil(stats(ig).BoundingBox(1)+stats(ig).BoundingBox(3)),size(bwi,2))) ...
                = bwi(max(floor(stats(ig).BoundingBox(2)),1):min(ceil(stats(ig).BoundingBox(2)+stats(ig).BoundingBox(4)),size(bwi,1)), max(floor(stats(ig).BoundingBox(1)),1):min(ceil(stats(ig).BoundingBox(1)+stats(ig).BoundingBox(3)),size(bwi,2)));
            bwp = bwperim(bw1reg);
            [r,c] = find(bwp==1,1);
            cellcontour = bwtraceboundary(bwi,[r c],'W',8,Inf,'counterclockwise');
        else 
            cellcontour =zeros(1,2);
             bw1reg = zeros(size(bwi));
        end
    end

    function bwo = segment(bwi)
        bwi = bwi./max(max(bwi));
        th = graythresh(bwi)+thresh;
        th = max(0,min(th,1));
        bwo = im2bw(bwi,th);       
        bwo = imclearborder(bwo);
        bwo = bwareaopen(bwo,5000);
        se = strel('disk',10);
        bwo = imclose(bwo,se);
        bwo = imfill(bwo,'holes');
        
    end

    function bwo = openi(Im)
        f1 = fspecial('Gaussian', openI, openI/3);
        f2 = fspecial('Gaussian', openI, openI/2);
        df = f1-2.*f2;
        bwo = conv2(double(Im), df, 'same');
        bwo = imcomplement(bwo);
    end

   
end