image thumbnail
from MATLAB codes for canopy image analysis by Lauri Korhonen
Codes for thresholding skyward-looking canopy images and calculating proportion of canopy pixels.

bin_cc.m
%
%	A MATLAB program for analysis of thresholded (=binary) canopy images
%
% 	Authors: Jaakko Heikkinen & Lauri Korhonen 2006-2009			
%	Citation: Korhonen, L. & Heikkinen, J. 2009. Automated analysis of in situ canopy images for the estimation of forest canopy cover.
%		Forest Science 55(4): 323-334.
%	Morphological analysis and rectangular AOV reduction for thresholded binary images  
%

%Constants

filetype='bmp';	%Which type of files are read
elem_size=10	%Diameter of the circular structuring element used in painting
true_focal=6.5;	%Focal length of the camera, mm
ccd_hor=6.16;	%Horizontal size of the ccd sensor
ccd_ver=4.62;	%Vertical size of the ccd sensor
hor_angle=90 	%Horizontal angle of view that should be used, degrees. Larger than cameras AOV -> whole image is analyzed.

files=dir(['*.', filetype]); 	%Read all files of the specified type

results=cell(length(files),3);	%Create array for the results
results{1,1}='Image';
results{1,2}='C_Closure'; 
results{1,3}='C_Cover'; 

disp('Analyzing image:');

%%%%%%%%%%%%%%% START LOOP %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

for h=1:length(files) %Loop goes through all images h

	%%%%%%%%%%%%%%%%%% READ IMAGE %%%%%%%%%%%%%%%%%%%%%%%%%%%%%

	fileid=files(h).name(1:length(files(h).name)-(length(filetype)+1));
	disp(fileid);
	
	binim=imread(files(h).name); 		%Read RGB image
	[rows,cols,channels]=size(binim);	%Read image dimensions

	%%%%%%%%%%%%% SCALE IMAGE TO NEW ANGLE OF VIEW %%%%%%%%%%%%%%%%

	%Calculate horizontal angle of view
	orig_angle=2*atan(ccd_hor/(2*true_focal))*180/pi;
	
	if orig_angle > hor_angle
		%Calculate the dimensions of the area to be analyzed
		new_cols=round(hor_angle*cols/orig_angle);
		new_rows=round(rows*new_cols/cols);
		
		cmin=round(cols/2-new_cols/2);
		cmax=round(cols/2+new_cols/2);
		rmin=round(rows/2-new_rows/2);
		rmax=round(rows/2+new_rows/2);
		
		binim=binim(rmin:rmax,cmin:cmax);	%Extract the area to be analyzed
		[rows,cols]=size(binim); 			%Update image dimensions
	
	end
	
	%%%%%%%%%%%%%%% CANOPY CLOSURE %%%%%%%%%%%%%%%%%%%%%%%%%%%%
	
	binim=binim/255;			%Scale 255 to 1 for calculations
	whites=sum(sum(binim));
	cclosure=(1-whites/(rows*cols))*100;

	%%%%%%%%%%%%%%%%%  CANOPY COVER %%%%%%%%%%%%%%%%%%%%%%%%%%%
	
	%Paint the crowns black using morphological opening and closing
	elem=strel('disk',elem_size);
	opeclo=imopen(binim,elem);	%Because objects are black, do the closing with opening operation
	opeclo=imclose(opeclo,elem);

	%Calculate pixel sum
	whites=sum(sum(opeclo));
	ccover=(1-whites/(rows*cols))*100;

	%%%%%%%%%%%%%% SAVE PAINTED IMAGE AND SPREADSHEET %%%%%%%%%%%%%%%%%%%%%

	%Save painted images as bmp
	opeclo=opeclo>0;	%2-bit conversion
	imwrite(opeclo, ['painted_',fileid,'_',num2str(hor_angle),'_deg.bmp'],'bmp')

	close all
	
	%Save results to array
	results{h+1,1}=files(h).name;
	results{h+1,2}=cclosure;
	results{h+1,3}=ccover;
	
end
results
xlswrite(['Binpainter',num2str(hor_angle),'.xls'],results);

clear all

Contact us at files@mathworks.com