Version 2

! Syntex of cGLOFLLS methods are changed. Please see the code and the example scripts. Image noise sensitivity analysis functions are added. Version and hardware will be checked when the processing option is 'gpu'. Downsampling(upsampling) process is improved. Graphical interface on the calibration process is improved. B16 image format is supported. Example scripts are moved into the main folder 'GLOFSFE'.
technolojin · Jun 25, 2019 · edbe681 · edbe681
1 parent f9bdd70
commit edbe681
Show file tree

Hide file tree

Showing 66 changed files with 2,668 additions and 722 deletions.
diff --git a/.gitignore b/.gitignore
@@ -1,5 +1,7 @@
-/data/*
+/GLOFSFE/Examples/DATA/*
 
 *.asv
 *.mat
 *.tif
+*.avi
+*.mp4
diff --git a/GLOFSFE/+SSA/horn_schunck.m b/GLOFSFE/+SSA/horn_schunck.m
@@ -1,6 +1,6 @@
 function [ Ux, Uy, error ] = horn_schunck( In1, In2, mask_1st, mask_2nd,...
     lambda, maxnum)
-%HORN_SCHUNCK horn-schunck optical flow estimator
+% HORN_SCHUNCK horn-schunck optical flow estimator
 % 
 % INPUT:
 %   In1, In2: image matrix 
@@ -13,17 +13,23 @@
 %   Ux, Uy: estimated optical flow
 %
 %
-%See also:
-% <a href="https://github.com/Tianshu-Liu/OpenOpticalFlow">
-%  OpenOpticalFlow</a>
-% <a href="http://doi.org/10.5334/jors.168"> Tianshu Liu, "OpenOpticalFlow:
-%  An Open Source Program for Extraction of Velocity Fields from Flow 
-%  Visualization Images"</a>
+% See also:
 %
+% * <https://github.com/Tianshu-Liu/OpenOpticalFlow 
+%    OpenOpticalFlow>
+%
+% * <http://doi.org/10.5334/jors.168
+%    Tianshu Liu, 
+%    "OpenOpticalFlow: An Open Source Program for Extraction of Velocity 
+%    Fields from Flow Visualization Images" 
+%    Journal of Open Research Software, 5(1), p.29. (2017)>
 %
-% Original work Copyright (c) 2017 Tianshu Liu
-% Released under the MIT license
 %
+% Original work 
+% Copyright (c) 2017 Tianshu Liu
+% Released under the MIT license
+% 
+% Modified:
 % Copyright (c) 2018 Taekjin Lee
 % Released under the MIT license
 % http://opensource.org/licenses/mit-license.php

diff --git a/GLOFSFE/+SSA/license.txt b/GLOFSFE/+SSA/license.txt
@@ -0,0 +1,22 @@
+The MIT License (MIT)
+
+Copyright (c) 2017 Tianshu Liu
+Copyright (c) 2018 Taekjin Lee
+
+Permission is hereby granted, free of charge, to any person obtaining a copy
+of this software and associated documentation files (the "Software"), to deal
+in the Software without restriction, including without limitation the rights
+to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
+copies of the Software, and to permit persons to whom the Software is
+furnished to do so, subject to the following conditions:
+
+The above copyright notice and this permission notice shall be included in all
+copies or substantial portions of the Software.
+
+THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
+OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
+SOFTWARE.
diff --git a/GLOFSFE/+SSA/liu_shen.m b/GLOFSFE/+SSA/liu_shen.m
@@ -14,12 +14,16 @@
 %   Ux, Uy: estimated optical flow
 %
 %
-%See also:
-% <a href="https://github.com/Tianshu-Liu/OpenOpticalFlow">
-%  OpenOpticalFlow</a>
-% <a href="http://doi.org/10.5334/jors.168"> Tianshu Liu, "OpenOpticalFlow:
-%  An Open Source Program for Extraction of Velocity Fields from Flow 
-%  Visualization Images"</a>
+% See also:
+%
+% * <https://github.com/Tianshu-Liu/OpenOpticalFlow 
+%    OpenOpticalFlow>
+%
+% * <http://doi.org/10.5334/jors.168
+%    Tianshu Liu, 
+%    "OpenOpticalFlow: An Open Source Program for Extraction of Velocity 
+%    Fields from Flow Visualization Images" 
+%    Journal of Open Research Software, 5(1), p.29. (2017)>
 % 
 % 
 % Original work:

diff --git a/GLOFSFE/+plot/imgQuantile.m b/GLOFSFE/+plot/imgQuantile.m
@@ -0,0 +1,7 @@
+function [ imgq ] = imgQuantile( img, Q)
+% 'Q' quantile of the given image
+temp=sort(nonzeros(img(:)));
+nt=size(temp,1)-1;
+imgq=temp(floor(nt.*Q)+1)';
+end
+
diff --git a/GLOFSFE/+plot/quiver_mod.m b/GLOFSFE/+plot/quiver_mod.m
@@ -7,6 +7,9 @@
 	offset = 1;
 end
 
+length = 1.5;
+
+
 [ny,nx]=size(U);
 [X,Y]=meshgrid(1:nx,1:ny);
 
@@ -33,16 +36,18 @@
 
 % Quantile of 'stdp' to the range
 stdp=0.95;
-temp=sort(nonzeros(I(:)));
-nt=size(temp,1);
-v_max=temp(floor(nt*stdp))/stdp;
+% temp=sort(nonzeros(I(:)));
+% nt=size(temp,1);
+% v_max=temp(floor(nt*stdp))/stdp;
+imgq = plot.imgQuantile(I, stdp);
+v_max=imgq/stdp;
 
 nk=size(C,1);
 Ic = round(I/v_max*(nk-1))+1;
 Ic(Ic>=nk)=nk;
 
-Unor=U./I*2*skp;
-Vnor=V./I*2*skp;
+Unor=U./I*length*skp;
+Vnor=V./I*length*skp;
 
 % mag=1.5*skp/v95;
 
@@ -51,6 +56,7 @@
 h=cell(nk,1);
 for k=1:nk
     idx=find(double(Ic(:)==k).*Mskip(:));
+    h{k}=quiver(X(idx),Y(idx),Unor(idx),Vnor(idx),0,'Color',[0,0,0],'LineWidth',2);
     h{k}=quiver(X(idx),Y(idx),Unor(idx),Vnor(idx),0,'Color',C(k,:),varargin{:});
 end
 set(gca,'Color','black','XLim',[1,nx],'YLim',[1,ny]);

diff --git a/GLOFSFE/+readB16/license.txt b/GLOFSFE/+readB16/license.txt
@@ -0,0 +1,24 @@
+Copyright (c) 2016, Carl Hall
+All rights reserved.
+
+Redistribution and use in source and binary forms, with or without
+modification, are permitted provided that the following conditions are
+met:
+
+    * Redistributions of source code must retain the above copyright
+      notice, this list of conditions and the following disclaimer.
+    * Redistributions in binary form must reproduce the above copyright
+      notice, this list of conditions and the following disclaimer in
+      the documentation and/or other materials provided with the distribution
+
+THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
+AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
+IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
+ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE
+LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
+CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
+SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
+INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
+CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
+ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
+POSSIBILITY OF SUCH DAMAGE.
diff --git a/GLOFSFE/+readB16/readB16.m b/GLOFSFE/+readB16/readB16.m
@@ -0,0 +1,109 @@
+% readB16.m
+function image = readB16(dirPath, filename)
+% readB16.m  Read PCO *.b16 image file
+%
+% image = readB16(filename)
+%       = readB16(dirPath, filename)
+%
+% where: image - 2D array of image
+%     filename - path and filename (string)
+%      dirPath - optional directory path (string)
+% 
+% Example usage: 
+%  For a file 'image.b16' located in the local directory 'files', the following
+%  will be successful
+%    image = readB16('files/image.b16')
+%    image = readB16('files/','image.b16')
+%
+%  and the following will cause failure with fd=-1
+%    image = readB16('files/image.B16')
+%    image = readB16('files','image.b16')
+%    image = readB16('files','image.B16')
+% 
+% Image is returned as a 2D array, where only 1 image is stored in each *.b16 file.
+%
+% Only b/w images have been tested, and importing color images should 
+% throw an error.
+%
+
+% Author: Carl Hall
+% Date: Spring 2016
+
+% Modification History
+%  March 2016 - Original Code
+%  June  2019 - Modified by Taekjin Lee
+
+%% b16 Header byte Structure
+% 
+% Bytes    Type      Description
+% 4        chars     "PCO-" 
+% 4        int32     File size in byte       
+% 4        int32     file size in byte 
+% 4        int32     header size + comment filed in byte
+% 4        int32     image width in pixel 
+% 4        int32     image height in pixel 
+% 4        int32     -1 (true), extended header follow
+% 4        int32     0 = black/with camera, 1 = color camer
+% 4        int32     black/white LUT-setting, minimum value 
+% 4        int32     black/white LUT-setting, maximum value 
+% 4        int32     black/white LUT-setting, 0 = linear, 1 = logarithmic 
+% 4        int32     red LUT-setting, minimum value 
+% 4        int32     red LUT-setting, maximum value 
+% 4        int32     green LUT-setting, minimum value 
+% 4        int32     green LUT-setting, maximum value 
+% 4        int32     blue LUT-setting, minimum value 
+% 4        int32     blue LUT-setting, maximum value 
+% 4        int32     color LUT-setting, 0 = linear, 1 = logarithmic 
+% ?        int32     internal use
+% ?        chars     Comment file in ASCII characters with variable length of 0...XX. 
+%                              The length of the comment filed must be documented in the “header length” field. 
+%          uint16    16 bit pixel data, starting at offset given by the 'header size' int32
+
+%% Start of Code
+
+% Parse optional input
+if nargin>1
+    filename = fullfile(dirPath,filename);
+else
+    filename = dirPath;
+end
+
+% Open the file
+fd = fopen(filename,'r');
+if(fd < 0)
+  error('Could not open file: %s',filename)
+end
+
+% Check that it is a PCO file
+filetype = fread(fd, 4, 'char');
+if(char(filetype') ~= 'PCO-')
+  error('Wrong filetype: %s',char(filetype'))
+end
+
+% Get the image dimensions:
+fileSize   = fread(fd, 1, 'int32');  % not used
+headLength = fread(fd, 1, 'int32');  % offset for image data
+imgWidth   = fread(fd, 1, 'int32');  %
+imgHeight  = fread(fd, 1, 'int32');  %
+
+
+% look into the extended header, and thow error if color image
+extHeader  = fread(fd, 1, 'int32');
+if(extHeader == -1)
+  colorMode  = fread(fd, 1, 'int32');
+  if(colorMode ~= 0)
+    error('Color image detected. Only b/w images have been tested with this function')
+  end
+end
+
+% Get the image
+fseek(fd, headLength, 'bof');
+image = fread(fd, [imgWidth,imgHeight], 'uint16');
+
+% rotate and flip image to suit user
+image = image';  
+
+% close the file
+fclose(fd);
+
+end
diff --git a/GLOFSFE/@cGLOFCase/cGLOFCase.m b/GLOFSFE/@cGLOFCase/cGLOFCase.m
@@ -15,7 +15,7 @@
 
         % data
         % calibration images
-        % 1:dark 2:bg 3:exc 4:scale 5:alpha
+        % 1:dark 2:bg 3:exc 4:scale 5:alpha(droplet)
         DirCal         % directories of calibration images folder, struct
 
         cal_fmt        % cal image format
@@ -90,20 +90,10 @@ function setFileMask(obj,file_mask)
         end
 
         function s = saveobj(obj)
-            s.Name = obj.Name;
-            s.Conditions = obj.Conditions;
-
-            s.DirCal=obj.DirCal;
-            s.cal_fmt=obj.cal_fmt;
-            s.max_image=obj.max_image;
-            s.FileMask=obj.FileMask;
-
-            s.scale_points=obj.scale_points;
-            s.scale_length=obj.scale_length;
-            s.flagScale=obj.flagScale;
-            s.oil_drops=obj.oil_drops;
-            s.v_drop=obj.v_drop;
-            s.flagDrops=obj.flagDrops;
+            f = properties(obj)';
+            for n=1:size(f(:),1)
+                s.(f{n})=obj.(f{n});
+            end
         end
         function save(obj,fname,varargin)
             s=obj.saveobj; %#ok<*NASGU>
@@ -116,21 +106,10 @@ function save(obj,fname,varargin)
         function obj = loadobj(s)
             if isstruct(s)
                 obj=cGLOFCase(s.Name);
-                obj.Name=s.Name;
-                obj.Conditions=s.Conditions;             
-
-                obj.DirCal=s.DirCal;
-                obj.cal_fmt=s.cal_fmt;
-                obj.max_image=s.max_image;
-                obj.FileMask=s.FileMask;
-
-                obj.scale_points=s.scale_points;
-                obj.scale_length=s.scale_length;
-                obj.flagScale=s.flagScale;
-                obj.oil_drops=s.oil_drops;
-                obj.v_drop=s.v_drop;
-                obj.flagDrops=s.flagDrops;
-
+                f=fieldnames(s)';
+                for n=1:size(f(:),1)
+                    obj.(f{n})=s.(f{n});
+                end
             else
                 obj=s;
             end

diff --git a/GLOFSFE/@cGLOFCase/setOilDrops.m b/GLOFSFE/@cGLOFCase/setOilDrops.m
@@ -1,26 +1,46 @@
-function setOilDrops(obj,oil_drops,v_drop)
+function setOilDrops(obj,v_drop,oil_drops)
+narginchk(1,3);
+
 if nargin==3
-    obj.oil_drops=oil_drops;
     obj.v_drop=v_drop;
-else
+    obj.oil_drops=oil_drops;
+elseif nargin<3
+    if nargin==2
+        obj.v_drop=v_drop;
+    else
+        volume_drop = input('please input the volume of an oil drop [micro liter]:');
+        obj.v_drop = volume_drop*10^-9; %[liter]
+    end
     img=cGLOFImageSet(obj.DirCal.alpha,obj.cal_fmt,obj.max_image);
     I_alpha=getIave(img);
-    imagesc(I_alpha,[0,1]);
-    n_drops = input('input number of oil drops:');
-    volume_drop = input('volume of an oil drop [micro liter]:');
-    volume_drop=volume_drop*10^-9; %[liter]
+    imgq=plot.imgQuantile( I_alpha, [0.05, 0.98]);
+    imagesc(I_alpha,imgq);
+    colormap('winter');
+    n_drops = input('please input the number of oil drops:');
     dp=zeros(n_drops,4);
     for i=1:n_drops
-        disp(['positon of dorp number ',int2str(i)]);
-        p = ginput(2);
+        disp(['positon of drop # ',int2str(i),'/',int2str(n_drops)]);
+        p=zeros(2);
+        hold on;
+        for n=1:2
+            p(n,:)=ginput(1);
+            scatter(p(n,1),p(n,2),'+r','SizeData',50,'LineWidth',1.5);
+        end
+        hold off;
         % Get the x and y corner coordinates as integers
         dp(i,1) = min(floor(p(1)), floor(p(2))); %xmin
         dp(i,2) = min(floor(p(3)), floor(p(4))); %ymin
         dp(i,3) = max(ceil(p(1)), ceil(p(2)));   %xmax
         dp(i,4) = max(ceil(p(3)), ceil(p(4)));   %ymax
+        hold on;
+        plot([dp(i,1),dp(i,1),dp(i,3),dp(i,3),dp(i,1)],...
+            [dp(i,2),dp(i,4),dp(i,4),dp(i,2),dp(i,2)],...
+            '-r','LineWidth',0.9);
+        text(dp(i,1)-10,dp(i,2),num2str(i,'%02u'),'BackgroundColor',[1,1,1],'Margin',0.1,'FontSize',8);
+        hold off;
     end
     obj.oil_drops=dp;
-    obj.v_drop=volume_drop;
+    pause(0.1);
 end
 obj.flagDrops=true;
 end