8
\$\begingroup\$

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.

ferada
11.4k25 silver badges65 bronze badges
asked Nov 30, 2015 at 21:22
\$\endgroup\$

1 Answer 1

5
\$\begingroup\$

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)
  • You're using spaces - Good!
    • (Most of the places, [pathstr,name,ext]=fileparts(filename);)
  • You're not using i and j 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.

answered Jul 12, 2016 at 14:34
\$\endgroup\$

Your Answer

Draft saved
Draft discarded

Sign up or log in

Sign up using Google
Sign up using Email and Password

Post as a guest

Required, but never shown

Post as a guest

Required, but never shown

By clicking "Post Your Answer", you agree to our terms of service and acknowledge you have read our privacy policy.

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.