%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% function process(image, foundtiles, theta, path) % out = process(image, foundtiles, theta) % image : grayscale image to be processed % foundtiles : tiles already located by user in first frame % theta : angle of rotation determined from first frame (ccw, degrees) % globals: OPTIONS - list of options, MASK - image to convolve with % % processed the grayscale image to try and find centers of tiles % method: % rotate image so that crystal is vertically oriented % convolve image with mask determined from first image % use tophat filter to emphsize small changes (option 4) % centers are determined by riongal maxima % nearby centers are combined global OPTIONS global MASK; [msizey, msizex] = size(MASK); yrange = floor(msizey/2); xrange = floor(msizex/2); % rotate image by determined angle ---------------------- % this is a very time consuming step rotim = imrotate(image, theta, 'bicubic'); figure(1) imagesc(rotim) colormap gray if ( ~isempty(find(OPTIONS == 12)) | ~isempty(find(OPTIONS == 11))) % allow user to change rotation angle for each frame -------- disp( 'create a vector points in the direcetion of the crystals orientation' ) disp( 'tiles should be oriented vertically') disp('create orientaton vector (2clicks) or keep old (center clicks)') [x,y,button] = ginput(2); if (button == 1) temp = atand((x(2)-x(1)) / (y(1)-y(2))); theta = theta + temp; rotim = imrotate(image, theta, 'bicubic'); figure(2) imagesc(rotim); end end [ymax,xmax] = size(rotim); save( path, 'rotim', '-mat', '-append'); save( path, 'theta','xmax','ymax', '-mat', '-append'); if ( isempty(find(OPTIONS == 11)) ) % rotate found tiles along with image temp = zeros(size(image)); for i=1:length(foundtiles) temp(foundtiles(i,2), foundtiles(i,1)) = 1; end temp = bwmorph( temp, 'dilate'); temp = imrotate(temp, theta); temp = bwmorph(temp , 'shrink', 2); [y,x] = find(temp); foundtiles = [x,y]; %convolve mask with image ----------------------------------- convim = conv2( double(rotim), MASK, 'same' ); % flatten black spots -------------------------------------- th = graythresh(convim); high = max(max(convim)); temp = (convim > high *th); convim = convim .* temp; % top hat -------------------------------------------------- se = strel( 'rectangle', [floor(msizey/2), floor(msizex/2)]); topim = imtophat(convim,se); % remove high peaks th = graythresh (topim); high = max(max(topim)); temp = (topim < th*high); topim = topim .* temp; % average results (remove tiny spots that are not really centers) h = fspecial('average', [yrange, xrange]); topim = filter2(h, topim); % average again to help combine centers later topim = filter2(h, topim); % mask out things that werent even in orig image temp = ones(size(image)); temp = imrotate(temp, theta); temp = bwmorph(temp, 'erode', xrange); topim = temp .* topim; % turn into a b/w image. hopefully with only areas round centers chosen--- if ( ~isempty(find(OPTIONS == 5)) ) % threshold image [i,j] = size(topim); thres = graythresh(topim( (i/2):(i/2+10*msixey), (j/2):(j/2+10*msixex)) ); high = max(max(topim)); bwim = (topim > thres*high); else % local max bwim = imregionalmax( topim ); end % mask out things not in orig image, and edge effects temp = bwmorph(temp, 'erode'); bwim = temp .* bwim; % remove vertical bars that result from top hat averaging bwim = bwmorph(bwim, 'shrink'); temp = bwmorph(bwim, 'clean'); bwim = (bwim & ~temp); % combine centers % this is a time consuming step. bwim = imdilate( bwim, ones(yrange, floor(msizex/3) )); bwim = bwmorph( bwim, 'thin', 7); % more? bwim = bwmorph( bwim, 'shrink', 10); temp = bwmorph(bwim, 'clean'); bwim = (bwim & ~temp); [y,x] = find(bwim); N = numel(x); else y = []; x = []; N = 0; end conn = zeros(N,6); save(path, 'x','y', 'N','conn','-mat', '-append'); % allow user to modify centeres if option is turned on if ( ~isempty(find(OPTIONS == 7)) ) modify( path ); end if ( ~isempty(find(OPTIONS == 9)) ) hairpin(path); else hair = zeros(size(x)); save(path, 'hair','-mat', '-append'); end if ~isempty(find(OPTIONS == 11)) modify( path ); % need to find x and y spacing (for connections) load('-mat', path, 'x','y','N'); ratio = 1.5; x = sqrt( N / ratio ); % 1.5 is assumption about rows vs columns NUM = ceil( 4 * x ); NUM2 = ceil( 4 * ratio * x ); % finding y distance ------------------------------------------------------ n = 1; for i = 1:N yi = y(i); for j = 1: (i-1) XY(n) = abs( yi - y(j) ); n = n+1; end end yhist = hist( XY, NUM ); % find max points of histogram high = zeros(5,1); j = 1; for i= 2:NUM if ( yhist(i-1) < yhist(i) & yhist(i) > yhist(i+1) ) high(j) = i; j = j + 1; if j == 6 break end end end % compute average distance yave = 0; for i = 1:4 yave = yave + high(i+1) - high(i); end yave = yave/4; yave = yave * imax / NUM % finding x distance ------------------------------------------------------ % if a point is within a certain y range, computer its distance to the % other points. % reuse = X vector range = yave / 3; n = 1; for i = 1:N xi = x(i); yi = y(i); for j = 1:i-1 if ( abs(yi - y(j)) < range ) XY(n) = abs( xi - x(j) ); n = n+1; end end end xhist = hist( XY(1:n) , NUM2 ); % find max points of histogram high = zeros(5,1); j = 1; for i= 2:NUM2 if ( xhist(i-1) < xhist(i) & xhist(i) > xhist(i+1) ) high(j) = i; j = j + 1; if j == 6 break end end end % compute average distance xave = 0; for i = 1:4 xave = xave + high(i+1) - high(i); end xave = xave/4; xave = xave * jmax / NUM2 save(path, 'xave', 'yave', '-mat', '-append'); end % %% MODIFY MASK % adds all tiles found to mask % empirically this doesnt perform all that well % [y,x] = find(bwim); % N =length (x); % tempMASK = zeros(size(MASK)); % for n = 1:N % tempMASK = tempMASK + rotim( (y(n)-yrange):(y(n)+yrange), (x(n)-xrange):(x(n)+xrange) ); % end % tempMASK = tempMASK / N; % tempMASK = tempMASK + fliplr(tempMASK) + flipud(tempMASK) + rot90(tempMASK,2); % tempMASK = tempMASK / 4; % % MASK = ( MASK + tempMASK ) / 2;