The code below takes a video as input (a drop onto a fabric surface) and I want it to return the area of the drop.
Currently the code works, but I'm no Matlab expert, and wrote it entirely through google/hack&slashing my way into a successful form. Unfortunately it's slow, and it's probably not written with good style.
I was hoping that anyone would be able to tell me how I can improve the code. I don't think that you'll need a video file or calibration file, but if it would help I can include it. I am also curious whether there is a better way to output the data to an XLS file than writing in the title/values.
I'm mainly concerned about style and efficiency.
%profile off
%profile clear
clear all;
set(0,'ShowHiddenHandles','on')
delete(get(0,'Children')) % Delete any currently open figures
%profile on
%Input user specified filename
[filename, pathname] = uigetfile('*.avi', 'Select a video file');
if isequal(filename,0)
disp('User selected Cancel')
else
disp(['User selected: ', fullfile(pathname, filename)]);
end
% [~,~,ext]=fileparts(filename);
dataFileName={filename};
nab=fullfile(pathname, filename);
[filename, pathname] = uigetfile({'*.xlsx'},'Select the calibration file'); % choose name and location of Excel calibration file
if isequal(filename,0)
disp('User selected Cancel');
else
disp(['User selected: ', fullfile(pathname, filename)]);
end
[pathstr,name,ext]=fileparts(filename);
calfile=fullfile(pathname, filename);
% read the calibration file specifically
interval = xlsread(calfile,'Calibration','B1');
repeats = xlsread(calfile,'Calibration','B2');
pixels_per_mm = xlsread(calfile,'Calibration','B8');
stainValue = xlsread(calfile,'Calibration','B10');
maxFrames = xlsread(calfile,'Calibration','B23');
maxstainNum = xlsread(calfile,'Calibration','B24');
densityBlood = xlsread(calfile,'Calibration','B41');
surftBlood = xlsread(calfile,'Calibration','B42');
densityAir = xlsread(calfile,'Calibration','B35') ;
visAir = xlsread(calfile,'Calibration','B36');
dsFrame = xlsread(calfile,'Calibration','B14') ; % frame for first dispensing of fluid
tsFrame = xlsread(calfile,'Calibration','B15') ; % frame for first image analysis tracking
% read the calibration file generally
[ndata, col1] = xlsread(calfile,'Calibration','A:A');
[col2, tt] = xlsread(calfile,'Calibration','B:B');
[ndata2, col3] = xlsread(calfile,'Calibration','C:C');
% Read the video file
spatterObj = VideoReader(nab);
nframes = spatterObj.NumberOfFrames; % determine the number of frames
if nframes > maxFrames % restrict the total number of frames to process to maxFrames
nframes = maxFrames;
end
I = read(spatterObj, 1);
%onlystains = zeros([size(I,1) size(I,2) 3 nframes], class(I)); % zero the array
onlystains = zeros([size(I,1) size(I,2) 3 nframes], 'like', I); % replace "class(I)" with I.
figure(100) % set up figure to plot
%hold on
x=0;
for k = tsFrame : nframes % process video file, frame by frame from the selected first frame
singleFrame = read(spatterObj, k);
% Convert to grayscale
I = rgb2gray(singleFrame);
I2 = imadjust(I);
% convert to binary image
%level = graythresh(I);
frame = im2bw(I2,stainValue);
% Invert image so stains are white objects on black background
frame=~frame;
% filter all objects with less than minNumpixels to clean background
minNumpixels = 500;
frame = bwareaopen(frame, minNumpixels);
% fill any holes in object(s)
frame = imfill(frame,'holes');
if k == tsFrame
displaystains = frame;
else
displaystains = displaystains + frame;
end
imshow(displaystains, 'InitialMagnification',100);
% locate stains in the kth frame
L = logical(frame);%convert array into logicals
if any(L(:))
props = regionprops(L, {'Centroid','FilledArea','Eccentricity','Orientation','EquivDiameter','MajorAxisLength','MinorAxisLength','Perimeter'});
end
numObj = numel(props); % number of stains in the frame
%nn(k)=numObj; % number of stains in kth frame
mj = 0;
parentStainMax=0;
% store stain data from each frame and plot in figure
%should preallocate size of the array with zeroes, to avoid
%copying etc
stainData = zeros(numObj,10);
for stain = 1 : numObj
if props(stain).Centroid(1) > 0 % can use whole image
stainData(stain,1) = props(stain).Centroid(1); % x-coord of stain
stainData(stain,2) = props(stain).Centroid(2); % y coord of stain
stainData(stain,3) = props(stain).FilledArea; % stain size (filled)
stainData(stain,4) = props(stain).Eccentricity; % ratio of distance between ellipse foci and major axis length, =0 for circle
stainData(stain,5) = props(stain).Orientation; % angle bewteen major axis of ellipse and x-axis
stainData(stain,6) = props(stain).EquivDiameter; % average diameter (non-filled)
stainData(stain,7) = props(stain).MajorAxisLength; % length of major axis
stainData(stain,8) = props(stain).MinorAxisLength; % length of minor axis
stainData(stain,9) = props(stain).Perimeter; % perimeter
stainData(stain,10) = k; % time - must be the last parameter in the list
end
if stainData(stain,3)> parentStainMax % find the parent stain (the largest object)
parentStainNo(k) = stain;
parentStainMax=stainData(stain,3);
parentStain(k)=stainData(stain,3)/(pixels_per_mm)^2; % store stain areas in mm^2
AR(k)=stainData(stain,8)/stainData(stain,7); % aspect ratio
peri(k)=stainData(stain,9)/(pixels_per_mm); % perimeter of stain in mm
circ (k)=peri(k)^2/(4*pi()*parentStain(k)); % compactness
K(k)=k; % frame number
end
end
end
X=parentStain;
%Ktrim=K(dsFrame:nframes);
figure(101)
h3 = plot(K,parentStain); % plot the change in area vs frame number
%hold off
% Invert Matrices
stainxCol = X.';
kCol = K.';
for n=1:nframes
if n < tsFrame
Time(n) = 0;
else
Time(n) = interval*(kCol(n)-dsFrame);
end
end
timeCol = Time.';
stainDiamCol = 2*sqrt(stainxCol/pi());
ARCol=AR.';
periCol=peri.';
circCol=circ.';
figure(101)
h3 = plot(Time,X); % plot the stain area vs time
% Write to excel spreadsheet
[filename, pathname] = uiputfile('*.xlsx','Save as'); % choose name and location of Excel file to save results in
stainresults = fullfile(pathname, filename);
message1 = 'Please wait while I write the data to Excel...';
h = msgbox(message1);
hExcel = actxserver('Excel.Application');
hExcel.SheetsInNewWorkbook = 1; % Limits number of sheets in about-to-be created workbook to 1
EW = hExcel.Workbooks.Add; % Starts a new workbook (not named yet)
EW.SaveAs(stainresults) % Name and close file
EW.Close
hExcel.Quit
Title1 = {'Calibration Data'};
Title2 = {'File Name: '};
xlswrite(stainresults,Title1,'Sheet1','A1')
xlswrite(stainresults,Title2,'Sheet1','B1')
xlswrite(stainresults,dataFileName,'Sheet1','C1')
xlswrite(stainresults,col1,'Sheet1','A3');
xlswrite(stainresults,col2,'Sheet1','B3');
xlswrite(stainresults,col3,'Sheet1','C3');
header3 = {'Frame', 'Time','Stain diameter','Stain Area', 'aspect ratio','perimeter','compactness' };
header4 = {'', '(sec)', '(mm)', '(mm^2)', '', '(mm)'};
sheetRef = 'stain ';
xlswrite(stainresults,header3,sheetRef,'A1'); % write the column headings
xlswrite(stainresults,header4,sheetRef,'A2');
comb = [kCol,timeCol,stainDiamCol,stainxCol,ARCol,periCol,circCol];
xlswrite(stainresults,comb,sheetRef,'A3'); % write the data
message2 = 'Finished writing data';
h = msgbox(message2);
%profile off
%profsave(profile('info'),'Profile_0.1_preallocating')
I've also included the pastebin link (http://pastebin.com/CCFcH893#) if you want syntax highlighting.
1 Answer 1
There's a lot going on in this script, but I've tried to cover all of the important parts.
Style and conventions:
I think your coding style is very good.
- You're using nice a descriptive variable names - Good!
- Your variable names are consitent (camelCase) - Good!
- (Except
pixels_per_mm
)
- (Except
- You're using spaces - Good!
- (Most of the places,
[pathstr,name,ext]=fileparts(filename);
)
- (Most of the places,
- You're not using
i
andj
as variable names in loops - Good! - You're using
.'
and not'
when transposing - Good!
Writing to Excel:
xlswrite
opens and closes the Excel server everytime you run the command. This can obiously be a huge bottleneck if you run the command in a loop, or many times in a row as you do. You have two alternatives, that are both quite good.
You are writing the individual headers first, followed by the individual columns with data. This can be done much simpler, either by first writing all headers, then all data like this:
headers = {'Calibration Data', 'Filename'}
allData = {col1, col2, col3}
xlswrite(stainresults, headers, sheetref, 'A1')
xlswrite(stainresults, allData, sheetref, 'A2')
or, put everything into a single cell array, and print it in one singe command:
headersAndData = {headers; allData};
xlswrite(stainresults, headersAndData, sheetref, 'A1')
The other alternative is to use xlswrite1
from the file exchange. It's a very good and versatile function that does everything xslwrite
does. The only difference is that xlswrite1
keeps Excel open until you're finished writing to it. This makes it a lot faster than xlswrite
when used in loops. This approach will be the simplest if you don't want to rewrite your code.
Reading from Excel:
You're reading one column at a time. Don't! Read all columns into cells in one go, and, if necessary split it into separate variables in Matlab afterwards:
[col1, col2, col3] = allData{:};
It's much better to keep it as one variable though and use indices when necessary. It's faster, and it's "the Matlab way".
Loops:
Loops such as this one should be rewritten:
for n=1:nframes
if n < tsFrame
Time(n) = 0;
else
Time(n) = interval*(kCol(n)-dsFrame);
end
end
First off: Pre-allocate memory for vectors that are used in loops, so in front of a loop like this you should always do: Time = zeros(nframes, 1);
. This is much faster. The Matlab editor even warns you about it. But, in this case, the loop isn't needed at all!
Time = interval.*(kCol - dsFrame);
Time(1:floor(tsFrame)) = 0; % or Time(1:tsFrame) = 0; if tsFrame is an integer
Extracting data from struct:
All the lines below can be simplified by converting props
into a cell
array:
stainData(stain,3) = props(stain).FilledArea;
stainData(stain,4) = props(stain).Eccentricity;
Simply do:
stainData = struct2cell(props(:))
This will create a cell array stainData
with numObj
columns, and rows corresponding to the fields in props
. You might need to tweak some other parts of the code to make it exactly the way you want it (since you have props.Centroid(1)
and props.Centroid(2)
, but that should be very simple.
Explore related questions
See similar questions with these tags.