% SEGMENTTIFS(VART,VARNM,VARZ,METHOD,SEGMCH,CELLDIAM,PLSTART,PLSTOP) 
% Load images from time start until time stop
% specified in a 1x2 vector VART. The wavelength is specified as strings in
% VARNM by using STRVCAT.
% Choices for segmentation method METHOD:
% METHOD = 1 : FINDEACHCELLSURF, using filling of original image to find minima
% METHOD = 2 : FINDSINGLECELLSURF, using thresholding to create minima
% METHOD = 3 : FINDEACHCELLSURFTH, using filling of threshold image to find
% METHOD = 4 : SEGMCT, finding cels in CT image
%
% SEGMCH is the number of the channel that shall be used for segmentation
% CELLDIAM is the expected cell diameter, varZ is the number of planes that
% are loaded. (This may not be changed!!)
% PLSTARt and PLSTOP are the numbering of the planes that shall be used for
% segmentation
function [] = segmcelltif(varT,varNM,method,segmCh,cellDiam,plStart,plStop)


% set to less than 1 is we want to resize image
resizePrm = 0.5;

test = 0;

%
% Find where the variation in the filename is
%
varZ = plStop-plStart + 1;
[ind,model] = findVariation(varT,varNM,varZ);


%
% Loop over all time stacks to load and segment
%


% number of char in time points
numCharT = ind(1,2) - ind(1,1) + 1;

% number of char in planes
numCharZ = size(int2str(varZ),2);

% number of channels
numCh = size(varNM,1);



for i = varT(1) : varT(2)

    disp('Loading images in stack');    

    % load all images in this stack
    for j = 1 : varZ
        next = 0;
        dpT = makeZero(numCharT,i);
        dpZ = makeZero(numCharZ,j);
        
        
        % load each channel
        for k = 1 : numCh
            nameIm = [model(1:ind(1,1)-1) dpT ...
                model(ind(1,2)+1 : ind(2,1)-1) varNM(k,:) ...
                model(ind(2,2)+1 : ind(3,1)-1) dpZ ...
                model(ind(3,2)+1 : end)]
            pause

            % try to read images
            try
                I(:,:,j,k) = double(imread(nameIm));
            catch
                next = 1;
                break;
            end;
            
            if next == 1
                break;
            end;            
                
           
        end;
        
    end;

    % use only seleected planes
    im = I(:,:,plStart:plStop,:);

    if resizePrm < 1
        for j = 1 : size(im,3)
            imTmp(:,:,j,1) = imresize(im(:,:,j,1),resizePrm,'bicubic');
            imTmp(:,:,j,2) = imresize(im(:,:,j,2),resizePrm,'bicubic');
        end;
        im = imTmp;clear imTmp;
    end;

    if test == 1
        show(im,1)
    end;

    % find cells
    if method == 1
        % finding each cell using filling
        [cellBW,wat] = findEachCellSurf(squeeze(im(:,:,:,segmCh)),cellDiam);
    elseif method == 2
        % the cellular area
        [cellBW,wat] = findSingleCellSurf(squeeze(im(:,:,:,segmCh)),cellDiam);
    elseif method == 3
        % finding each ccell using adaptive thresholding and filling!
        [cellBW,wat] = segmsurfwat(squeeze(im(:,:,:,segmCh)),cellDiam);
    elseif method == 4
        % CTG staining
        [cellBW,wat] = segmct(squeeze(im(:,:,:,segmCh)),cellDiam);
    else
        error('Wrong option');
    end;
       

    % save the segmentation
    name = ['im' int2str(i) 'SegmAut'];
    disp(sprintf('Saving %s',name))
    ver = version;
    if strcmp(ver(1:3),'6.5')
        save(name,'im','cellBW','wat'); 
    else
        save(name,'im','cellBW','wat','-V6'); 
    end;            

    if test == 1
        show(wat,2)   
        show(cellBW,3)
    
        pause
    end;

end;


%-------------------------------------------------

function [ind,model] = findVariation(varT,varNM,varZ);

%
% Find the indices for the variation
%
all = dir('*tif');
model = all(1).name;
mName = ones(1,size(model,2));

for j = 1 : size(all,1)
    a = all(1).name;
    b = all(j).name;
    for i = 1 : size(a,2)            
        m = strmatch(a(i),b(i));
        if  mName(i) == 1 & isempty(m)
            mName(i) = 0;
        end;                
    end;        
  
end;


% find the indices for start and stop for variation
ind = [];teller = 1;
for j = 1 : size(mName,2)-1
   if mName(j) == 1 & mName(j+1) == 0
       ind(teller,1) = j+1;
   elseif mName(j) == 0 & mName(j+1) == 1
       ind(teller,2) = j;
       teller = teller + 1;
   end;       
end;


% to account for equal last and second last can be equal, search for the NM indices
for i = 1 : size(varNM,1)
    indVarNM = findstr(varNM(i,:),model);
    if ~isempty(indVarNM)
        ind(2,1) = indVarNM;
        ind(2,2) = indVarNM + 2;
    end;
end;


%----------------------------------------

function [num] = makeZero(p,i)


sizeI = size(int2str(i),2);
zero = num2str(zeros(1,p - sizeI));
ind = isspace(zero);
zero(find(ind == 1)) = [];
num = [zero int2str(i)];





