Matlab Computer Vision

Question

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.

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)
• 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.