基于PCA 和迭代 Canny Edge皮肤病变分割算法研究附Matlab代码.zip

%% Quicklook script for mole detection in cellphone images clear all % Parameters out_on=0; save_on=1; % File tree walker % fname='internet/atypical_naevi.jpg'; % fname='Stanford/ssm/SSM.png'; % fname='Stanford/ssm/Image_01.jpg'; % fname='Stanford/ssm/Image_06.jpg'; % fname='Stanford/normal/CS_1-26-05 B.tif'; %fname='R:\Project\New folder (2)\benign1.bmp'; % fname='Stanford/normal-fp/CS_10-24-01.tif'; global im im2 img [path,user_cance]=imgetfile(); if user_cance msgbox(sprintf('Error'),'Error','Error'); return end im=imread(path); im=im2double(im); %converts to double img=im; %for backup process :) %axes(handles.axes3); imshow(im); % Load image %img=imread(fname); %img=imresize(img,512/size(img,1)); % Convert to grayscale imgbw=rgb2gray(img); % Display image if out_on figure(1) %subplot(2,2,1); imagesc(imgbw) colormap 'jet' title(fname,'interpreter','none') end %% Analyze % Binarize using Otsu's method img_bn=~(im2bw(imgbw,graythresh(imgbw))); % Skin lesions are darker if out_on figure(2) imshow(img_bn) title('Binarized image') end %% Binarize using locally adaptive thresholds % % % Parameters % tl_sz=64; % Tile Size % lv_tsh=0; % Local variance threshold % gl_tsh=graythresh(imgbw); % Uniform regions threshold (global threshold) % % % Initialize variables % nu_rgns=ones(size(imgbw,1)/tl_sz,size(imgbw,2)/tl_sz); % Non-uniform mask % tsh=zeros(size(imgbw,1)/tl_sz,size(imgbw,2)/tl_sz); % Local thresholds % imgbw=double(imgbw); % % % Apply local thresholding % for yy=1:size(imgbw,1)/tl_sz % for xx=1:size(imgbw,2)/tl_sz % % Extract a tile from the image % y_bnds=(yy-1)*tl_sz+1:yy*tl_sz; % x_bnds=(xx-1)*tl_sz+1:xx*tl_sz; % tl=imgbw(y_bnds,x_bnds); % % % Compute local variance % l_var=var(double(tl(:))); % % % Calculate threshold based on local variance % if l_var<lv_tsh % % Uniform region % nu_rgns(yy,xx)=0; % tsh(yy,xx)=gl_tsh; % else % % Non-uniform region % tsh(yy,xx)=graythresh(tl); % end % end % end % % % Binarize image based on local thresholding map % tsh2=imresize(tsh,[size(imgbw,1) size(imgbw,2)],'bilinear'); % nu_rgns2=imresize(nu_rgns,[size(imgbw,1) size(imgbw,2)],'nearest'); % img_b=imgbw>tsh2; % img_b(~nu_rgns2)=imgbw(~nu_rgns2)>gl_tsh; %% Small region removal img_bns=img_bn; %Parameters p_sz_thr=500; % Positive region threshold n_sz_thr=1500; % Negative region threshold % Small region removal on positive image img_lbl=bwlabel(img_bns,4); bins=1:max(img_lbl(:)); a=histc(img_lbl(:),bins); blist=bins(a<p_sz_thr); img_bns(ismember(img_lbl,blist))=0; % img_lbl=bwlabel(img_bns,4); % for rgn=1:size(unique(img_lbl),1) % idx=find(img_lbl==rgn); % if size(idx,1) < p_sz_thr % img_bns(idx)=0; % end % end % Small region removal on negative image img_lbl=bwlabel(~img_bns,4); bins=1:max(img_lbl(:)); a=histc(img_lbl(:),bins); blist=bins(a<n_sz_thr); img_bns(ismember(img_lbl,blist))=1; % img_lbl=bwlabel(~img_bns,4); % for rgn=1:size(unique(img_lbl),1) % idx=find(img_lbl==rgn); % if size(idx,1) < n_sz_thr % img_bns(idx)=1; % end % end if out_on figure(3) imshow(img_bns) title('Small region removal') end %% Identify primary region of interest % Look for the connected component closest to center of image img_x0=size(imgbw,2)/2; img_y0=size(imgbw,1)/2; img_lbl=bwlabel(img_bns,4); stats=regionprops(img_lbl,'Centroid'); mindist=inf; for rgn=1:numel(stats) dist=sqrt(((stats(rgn).Centroid(1)-img_x0)^2) + ... ((stats(rgn).Centroid(2)-img_y0)^2)); if dist<mindist mindist=dist; minrgn=rgn; end end img_pr=img_bns; img_pr(img_lbl~=minrgn)=0; if out_on figure(4) imshow(img_pr) title('Primary Region of Interest') end %% Boundary detection img_ed=edge(img_pr); [edgs_x edgs_y]=ind2sub(size(img_ed),find(img_ed)); if out_on figure(5) imshow(img_ed) title('Edge Detection') end %% Analyze for border irregularity % % Move origin to region centroid % cn_edgs_x=edgs_x-stats(minrgn).Centroid(2); % cn_edgs_y=edgs_y-stats(minrgn).Centroid(1); % % figure(7) % clf % hold on % plot(cn_edgs_x,cn_edgs_y,'.') % hold off % axis equal % % % Convert to polar coordinates % [th r]=cart2pol(cn_edgs_x,cn_edgs_y); % [th idx]=sort(th); % r=r(idx); % % figure(8) % plot(th,r) % % % Take power spectral density % h=abs(fft(r)).^2; % figure(9) % semilogy(h) %% Analyze for symmetry % Convert to x,y list idx=find(img_pr); [rgn_x rgn_y]=ind2sub(size(img_pr),idx); % Calculate symmetric indices rgn_x_rot=round(2*stats(minrgn).Centroid(2)-rgn_x); rgn_y_rot=round(2*stats(minrgn).Centroid(1)-rgn_y); idx_rot=sub2ind(size(img_pr),rgn_x_rot,rgn_y_rot); pt=1; if out_on figure(7) imshow(img) hold on plot(rgn_y(pt),rgn_x(pt),'.') plot(rgn_y_rot,rgn_x_rot,'r.') plot(stats(minrgn).Centroid(1),stats(minrgn).Centroid(2),'g.') hold off end % Calculate symmetry MSE A=mean((img_bns(idx)-img_bns(idx_rot))/2); %% Calculate image intensity variance within region C=std(double(imgbw(img_pr))); %% Analyze border strength r=1.0;g=1.0;b=0.0; D=sqrt((r-0.1^2)+(g-0.1^2)+(b-0.1^2)); % Calculate gradient magnitude sigma=2; G=fspecial('gaussian',30,sigma); [dx dy]=gradient(G); img_dx=imfilter(double(imgbw),dx,'symmetric'); img_dy=imfilter(double(imgbw),dy,'symmetric'); img_grad=sqrt(img_dx.*img_dx + img_dy.*img_dy); idx_ed=find(img_ed); B=mean(img_grad(idx_ed)); if out_on figure(8) imagesc(img_grad) colormap 'jet' hold on plot(edgs_y,edgs_x,'r.','MarkerSize',2) hold off end %% Annotate Original Image f=figure(6); imshow(img) hold on % plot(stats(minrgn).Centroid(1),stats(minrgn).Centroid(2),'rx', ... % 'MarkerSize',20,'LineWidth',3) plot(edgs_y,edgs_x,'r.') hold off title('Annotated Image') an_array={['Asymmetry MSE: ' num2str(A)],... ['Border Strength: ' num2str(B)],... ['Color Standard Deviation: ' num2str(C)],... ['Diameter or Differential Structure',num2str(D)]}; text(10,10,an_array,'VerticalAlignment','top','FontSize',20) TDS=A*1.3+B*0.1+C*0.5+D*0.5; if TDS>5.45 msgbox(sprintf('CANCER DETECTED'),'Complete','Complete'); %return; else msgbox(sprintf('NORMAL SKIN'),'',''); end; %%if save_on %% [~, name, ~]=fileparts(fname); %% print(f,'-djpeg','-r50',[name '.jpg']) %%end