function [cmcontrol,controlim] = controlroi(imsegm,labelim,nuclei,rim,prm)

msg = ['This is ' upper(mfilename) ' finding control regions'];
disp(msg);

tic
% find possible did cells (NB lower threshold than earlier)
bw = imsegm > prm.detect.leveldid;

% remove noise objects
bw = bwareaopen(bw,10);

% label and find CM
[faser,L] = bwlabeln(bw);

msg = ['Finding ' int2str(L) ' potential DiD objects'];
disp(msg);

cm = zeros(L,2);
for i = 1 : L
    % this region
    reg = faser == i;
        
    clear c;
    [c(:,1),c(:,2)] = find(reg);

    % center of mass
    cm(i,:) = mean(c,1);
end;
clear faser;
cm = round(cm);
toc;

% distance around cells
dist = prm.detect.smallestdistcells;

% save test
% pause

% label image around the DiD cells not to take these regions
dim = size(labelim);
controlim = zeros(dim(1),dim(2),'uint8');
for i = 1 : L
    x = max(1,cm(i,1) - dist):min(dim(1),cm(i,1)+dist);
    y = max(1,cm(i,2) - dist):min(dim(2),cm(i,2)+dist);
    controlim(x,y) = 1;    
end;
% show(controlim,1)

% % all possible positions to use as control image
% clear c;
% [c(:,1),c(:,2),indcontrolim] = find(controlim == 0);

% % order them randonly
% n = size(c,1);
% v = randperm(n);
% c = c(v,:);

% rim coordinates
[rim.c(:,1),rim.c(:,2)] = find(rim.bw);
rim = rmfield(rim,'bw');

msg = ['Stopping after finding ' int2str(prm.detect.ncontrol) ' control ROIs'];
disp(msg);

% count number of wells
[faser,nwells] = bwlabeln(labelim > 0);

% save test2

ncontroltotal = 0;
for j = 1 : nwells

    
    msg = ['Well number ' int2str(j)];
    disp(msg);
    
    % stay within this well
    reghere = faser == j;
    % and stay within the possible control pixels
    reghere = reghere .* (controlim == 0);

    clear cwell;
    [cwell(:,1),cwell(:,2)] = find(reghere);    
    ntot = size(cwell,1);
    % randomize order
    v = randperm(ntot);
    cwell = cwell(v,:);
    clear reghere;
    
    ncontrolwell = 0;
    cmcontrol = zeros(prm.detect.ncontrol,2);
    for i = 1 : ntot
        cmhere = cwell(i,:);
        x = max(1,cmhere(1)-dist):min(dim(1),cmhere(1)+dist);
        y = max(1,cmhere(2)-dist):min(dim(2),cmhere(2)+dist);
        reg = controlim(x,y);

        % is there an overlap to a DiD cell image?
        if ~isempty(find(reg,1))
            msg = ['Overlap to DID or earlier object'];
            disp(msg);        
            continue;
        end;

        % is it dense enough of cells?     
        imhere = nuclei.bw(x,y);
        [a,n] = bwlabeln(imhere);
        if n < prm.detect.minnumneigh
            msg = ['Too sparse cells'];
            disp(msg);
            continue;
        end;

        % is it far enough from rim?
        d = rim.c - repmat(cmhere.*prm.lowres.h,rim.nc,1);    
        d = d.^2;
        d = sum(d,2);
        d = sqrt(d);
        d = min(d);
        if d < prm.detect.smallestdistrim
            msg = ['Too close too rim'];
            disp(msg);
            continue;
        end;

        % label this region as "taken"
        controlim(x,y) = 2;

        % and store the coordinates
        ncontrolwell = ncontrolwell + 1;
        ncontroltotal = ncontroltotal + 1;
        cmcontrol(ncontroltotal,:) = cmhere;
        msg = ['Finding control ROI ' int2str(ncontrolwell)];
        disp(msg);

        if i == ntot
            msg = ['Could not find enough control regions'];
            disp(msg);
            break;
        end;
        
        if ncontrolwell == prm.detect.ncontrol
            msg = ['Reaching maximum number of control ROIS'];
            disp(msg);
            break;
        end;
        
    end;
end;




