%% the saved file "*_res_nucl.mat" contains the BW_golgi
%% and cell centroid. In this one, the coordinates are stored in (Y,X) WITH the 0 in top-left
%% meaning that BW_nucl(nucl_centroid_imgRef(1,i), nucl_centroid_imgRef(2,i), i) = 1 [counterintuitive]
function segment_nucleus(send_PathName,acquisition_name,BW_golgi,golgi_centroid_imgRef,BW_cell,cell_centroid_imgRef)
%% 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 = '';
if nargin >=1 
    PathName= send_PathName;
else
 PathName = '';
end

frame = 1;
Nf = 150;
Img = [];
BW_nucl = [];nucl_centroid_imgRef =[];
nucl_centroid_orthoRefX = [];
nucl_centroid_orthoRefY = [];
thresh_img = [];
thresh = 0.5;
trX = []; trY = [];
scale_param_pos = 1/0.2167;

%% 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','Proceed data','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',PathName,...
    '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',0.0,'Max',1.0, 'value',thresh,...
    'Position',[50,670,100,30],...
    'SliderStep',[0.4/200 0.1],'Callback',@hthresh_Callback);

hsegfow = uicontrol('Style','pushbutton','String', 'Seg >>',...
    'Position',[160,670,70,30],'FontSize',12,...
    'Callback',@hsegfow_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,...
    hframe,hsegfow,hNextSegbutton,...
    hthresh],...
    'Units','normalized');
%% Final settings
% Assign the GUI a name to appear in the window title.
set(segfigure,'Name','Segment Nucleus v 1.0')
% Move the GUI to the center of the screen.
movegui(segfigure,'center')
% Make the GUI visible.
set(segfigure,'Visible','on')


if nargin>=1  
    tiff_files = dir([PathName '\' acquisition_name '*.tif']); % look for the tiff files
    for i = 1:length(tiff_files)
        if( ~isempty(strfind(tiff_files(i).name,'_405.tif') ) )
            FileName = tiff_files(i).name;
        end
    end
    loadData();
end

%% Callbacks
    function keypress(~,~)
        % Not used...
    end
    function click(~,~)
        % Not used...
    end
    function hquitbutton_Callback(~,~) % ajouter le lancement de l'autre GUI, enregistrer les contours de chaque cellule
        
        isOK_closing = checkSegmentationForAllImages();
       
        if( isOK_closing == 1)
            computeNucleusCentroids();
            
            close(segfigure); 
            
            save(fullfile(PathName,[acquisition_name 'res_nucleus.mat']),'nucl_centroid_imgRef','BW_nucl');
        end
        
    end

    function isOK_closing = checkSegmentationForAllImages()

        isOK_closing = 1;

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


    function computeNucleusCentroids(~,~)
        h = waitbar(0,'Saving informations...');
        nucl_centroid_imgRef = zeros(2,Nf);
        for framenumber = 1:Nf
            waitbar(framenumber/Nf,h);
            [BW_nucl(:,:,framenumber),~] = computeCentroidOneImage(framenumber); 
        end
        delete(h);
    end

    function hNextSeg_Callback(~,~) % ajouter le lancement de l'autre GUI, enregistrer les contours de chaque cellule
        
        isOK_closing = checkSegmentationForAllImages();
        
        if( isOK_closing == 1)
            computeNucleusCentroids();
            
            save(fullfile(PathName,[acquisition_name 'res_nucleus.mat']),'nucl_centroid_imgRef','BW_nucl');
            close(segfigure); 
            visualize_centroids_trajectories(PathName,acquisition_name,BW_nucl,nucl_centroid_imgRef,BW_golgi,golgi_centroid_imgRef,BW_cell,cell_centroid_imgRef); %PathName_send,acquisition_name,BW_nucl,nucl_centroid_imgRef,BW_golgi,golgi_centroid_imgRef,BW_cell,cell_centroid_imgRef
        end 
        
    end

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

            cell_file = dir([PathName acquisition_name  '*res_cell.mat']);
            golgi_file = dir([PathName acquisition_name  '*res_golgi.mat']);
            if( isempty(cell_file) || isempty(golgi_file)  )
                msgbox('You should first segment the cell and golgi of your data set','Incorrect file','error');
            else
                data_cell = load(fullfile(cell_file.folder,cell_file.name));
                data_golgi = load(fullfile(cell_file.folder,golgi_file.name));
                checkCellGolgiMatFormat(data_cell,data_golgi);
                loadData();
                showI(Img(:,:,1),BW_nucl(:,:,1));
            end
        end
        
    end


    function checkCellGolgiMatFormat(data_cell,data_golgi) 
        % si telle varialbe existe la convertir
        if ( isfield(data_cell,'cell_centroid') ) % old version, put the result in "cell_centroid_imgRef"
            cell_centroid_imgRef = data_cell.cell_centroid;
        else
            cell_centroid_imgRef = data_cell.cell_centroid_imgRef;
        end
        
        if ( isfield(data_golgi,'golgi_centroid') ) % old version, put the result in "cell_centroid_imgRef"
            golgi_centroid_imgRef = data_golgi.golgi_centroid;
        else
            golgi_centroid_imgRef = data_golgi.golgi_centroid_imgRef;
        end
    end


    function [bw1reg,nuclcontour] = computeCentroidOneImage(framenumber)
        se = strel('disk',20);
        bwi = imclose(BW_nucl(:,:,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 probable" is chosen
                d=zeros(ss,1);
                for j=1:ss
                    %uses the previous golgi centroid except for the first
                    %image or if the previous nucleus was not segmented
                    if( framenumber == 1 || nucl_centroid_orthoRefX(framenumber-1) == 0 )
                        d(j)=(stats(j).Centroid(1)-cell_centroid_imgRef(2,framenumber)).^2+(stats(j).Centroid(2)-cell_centroid_imgRef(1,framenumber)).^2; % centroid in raw, column center_cell in x,y
                    else
                        centroidX_orthoRef = stats(j).Centroid(1)+trX(framenumber);
                        centroidY_orthoRef = sI(1)-stats(j).Centroid(2)+trY(framenumber);
                        d(j)=(centroidX_orthoRef-nucl_centroid_orthoRefX(framenumber-1)).^2+(centroidY_orthoRef-nucl_centroid_orthoRefY(framenumber-1)).^2;
                    end
                end
                [~,ig] = min(d);
            else
                ig=1;
            end
            nucl_centroid_imgRef(:,framenumber) = [stats(ig).Centroid(2),stats(ig).Centroid(1)];

            nucl_centroid_orthoRefX(framenumber) = nucl_centroid_imgRef(2,framenumber)+trX(framenumber);
            nucl_centroid_orthoRefY(framenumber) = sI(1)-nucl_centroid_imgRef(1,framenumber)+trY(framenumber);
            
            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);
            nuclcontour = bwtraceboundary(bwi,[r c],'W',8,Inf,'counterclockwise');
        else 
            nuclcontour =zeros(1,2);
        end
    end


    function loadData()
        info = imfinfo([PathName,FileName]);
        Nf = length(info);
        Wd = info(1).Width;
        Hg = info(1).Height;
        sI = [Hg,Wd];
        h = waitbar(0,'Please wait...');
        Img = ReadTiffStack([PathName,FileName]);

        BW_nucl = zeros(size(Img,1),size(Img,2),size(Img,3));
        thresh_img = thresh*ones(Nf,1);
        for ind_img=1:Nf
            waitbar(ind_img/Nf,h);
            
            I = Img(:,:,ind_img);
            
            I = imadjust(I); % increases contrast
            %% first mask: cell border
            bw = im2bw(I,thresh);

            %bw = bwareaopen(bw,50);
            % to obtain a better shape of the object
            se = strel('disk',20);
            bw = imclose(bw,se);
            bw = imfill(bw,'holes');
            
            BW_nucl(:,:,ind_img) = bw;
        end

        delete(h);
        
        set(hframe,'max',Nf,'SliderStep',[1/(Nf-1) 0.2])
        set(hmessages,'String',['Loaded. Frame ',num2str(frame)]);
        compute_shifts_on_Image();
        
        nucl_centroid_imgRef = zeros(2,Nf);
        nucl_centroid_orthoRefX = zeros(1,Nf);
        nucl_centroid_orthoRefY = zeros(1,Nf);
    end


    function hframe_Callback(~,~)
        frame = get(hframe,'Value');
        frame = round(frame);
        set(hframe,'Value',frame);
        
        showI(Img(:,:,frame),BW_nucl(:,:,frame));
        set(hmessages,'String',['Frame ',num2str(frame)]);
        set(hthresh,'Value',thresh_img(frame,1));
    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_nucl(:,:,frame));
        set(hmessages,'String',['Frame ',num2str(frame)]);
        set(hthresh,'Value',thresh_img(frame,1));
        
    end

    function hthresh_Callback(~,~)
        thresh = get(hthresh,'Value');
        thresh_img(frame,1) = thresh;
        BW_nucl(:,:,frame) = segment(Img(:,:,frame));
        showbw(BW_nucl(:,:,frame));
        showI(Img(:,:,frame),BW_nucl(:,:,frame));
    end

    function hsegfow_Callback(~,~)
        h = waitbar(0,'Please wait...');
        for ind_f = frame:Nf
            waitbar(ind_f/Nf,h)
            thresh_img(ind_f,1) = thresh;
            BW_nucl(:,:,ind_f) = segment(Img(:,:,ind_f));
        end
        delete(h);
    end

%% Functions
    function showI(image,bwi) % ajouter la comparaison au centroid 
        axes(ha);
        imshow(image,[]); hold on;
        
        [~,nuclcontour] = computeCentroidOneImage(frame);
        plot(nuclcontour(:,2),nuclcontour(:,1),'r','LineWidth',0.5);

        hold off
        axes(hbw);
        imshow(bwi,[]);
    end

    function showbw(bwimage)
        axes(hbw);
        imshow(bwimage,[]);
    end

    function bwo = segment(I)
        tmp = double(I)/double(max(I(:)));
        bw= tmp>thresh;
        %bw = bwareaopen(bw,50);
        % to obtain a better shape of the object
        se = strel('disk',20);
        bw = imclose(bw,se);
        bwo = imfill(bw,'holes');

    end


    function compute_shifts_on_Image()
         pos_file = dir([PathName '\' acquisition_name '*Positions.mat']);
    
         if ( ~isempty(pos_file) )
             load(fullfile(PathName,pos_file.name));
         else
             warning = 'No positions file found, impossible to compute the translations';
         end
         
        X = p(1,:)-p(1,1);
        Y = p(2,:)-p(2,1);
        
        trX = X*scale_param_pos; % carefull there seems to be more positions than the number of images
        trY = Y*scale_param_pos;
         
    end

   
end