Skip to content

francesco changes in cross_correlation_rect #19

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
wants to merge 3 commits into
base: master
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
7 changes: 5 additions & 2 deletions cross_correlate_rect.m
Original file line number Diff line number Diff line change
@@ -1,4 +1,4 @@
function [c] = cross_correlate_rect(a2,b2,NfftHeight,NfftWidth)
function [c,c2] = cross_correlate_rect(a2,b2,NfftHeight,NfftWidth)

% temprorary solution
a2 = a2 - mean2(a2);
Expand All @@ -13,8 +13,11 @@

% another option not implemented yet, but b2 shall be larger than a2
% see >> help normxcorr2
a3=a2(2:end-1,2:end-1); %taglio bordi di a2
b2 = b2(end:-1:1,end:-1:1); %raddrizzo l'immagine
c2 = normxcorr2(a3,b2); %calcolo correlazione normalizzata
% c = normxcorr2(b2,a2);

c(c<0) = 0;

c2(c2<0) = 0;
return
4 changes: 3 additions & 1 deletion fill_holes.m
Original file line number Diff line number Diff line change
Expand Up @@ -7,7 +7,9 @@
%

vector(abs(vector)==0) = NaN;
vector = inpaint_nans(real(vector)) + i*inpaint_nans(imag(vector));
%vector = inpaint_nans(real(vector)) + i*inpaint_nans(imag(vector));
%%SOSTITUZIONE MIA
vector = inpaint_nans(vector,3);

% [indx,indy] = find(~abs(vector));
%
Expand Down
Binary file added lastpath.mat
Binary file not shown.
48 changes: 30 additions & 18 deletions openpivgui.m
Original file line number Diff line number Diff line change
Expand Up @@ -472,7 +472,6 @@ function pushbutton_start_Callback(hObject, eventdata, handles)

switch handles.filesType
case{'sequence'}

for fileind = 1:handles.amount-jump % main loop, for whole file list
% while (get(handles.pushbutton_start,'UserData') ==1)

Expand All @@ -493,9 +492,9 @@ function pushbutton_start_Callback(hObject, eventdata, handles)
% Prepare the results storage;
numcols = floor((horSize-ittWidth)/ovlapHor+1);
numrows = floor((verSize-ittHeight)/ovlapVer+1);
res = zeros(numcols*numrows,5);
res = zeros(numcols*numrows,6);
resind = 0;

%%%%
a2 = zeros(ittHeight,ittWidth);
b2 = zeros(ittHeight,ittWidth);
NfftWidth = 2*ittWidth;
Expand All @@ -516,7 +515,9 @@ function pushbutton_start_Callback(hObject, eventdata, handles)
% a2 = prepfun(a2);
% b2 = prepfun(b2);

c = cross_correlate_rect(a2,b2,NfftHeight,NfftWidth);
[c,c2] = cross_correlate_rect(a2,b2,NfftHeight,NfftWidth);

corr=max(max(c2));
% c = cross_correlate_rect(a2,b2,Nfftx,Nffty);
if ~any(c(:)), % completely "black"
u = 0;
Expand All @@ -525,7 +526,7 @@ function pushbutton_start_Callback(hObject, eventdata, handles)
x = origin(1) + k + ittWidth/2 - 1;
resind = resind + 1;
s2n = 0;
res(resind,:) = [x y u v s2n];
res(resind,:) = [x y u v corr s2n];
continue
end

Expand All @@ -540,7 +541,7 @@ function pushbutton_start_Callback(hObject, eventdata, handles)
x = origin(1) + k + ittWidth/2-1;

resind = resind + 1;
res(resind,:) = [x y u v s2n];
res(resind,:) = [x y u v corr s2n];
% quiver(x+cropvec(1),y+cropvec(2),u,v,'y');
if u ~= 0 || v ~= 0
% quiver(x,y,u,v,5,'y','Linewidth',1);
Expand All @@ -562,17 +563,22 @@ function pushbutton_start_Callback(hObject, eventdata, handles)
% Reshape U and V matrices in two-dimensional grid and produce
% velocity vector in U + i*V form (real and imaginary parts):

u = reshape(res(:,3), numrows,numcols);
v = reshape(res(:,4), numrows,numcols);
vector = u + sqrt(-1)*v;
u = reshape(res(:,3), numcols,numrows); %addition
v = reshape(res(:,4), numcols,numrows); %addition
u=u'; %addition
v=v'; %addition
vector=u + sqrt(-1)*v; %addition
% u = reshape(res(:,3), numrows,numcols);
% v = reshape(res(:,4), numrows,numcols);
% vector = u + sqrt(-1)*v;

% Remove outlayers - GLOBAL FILTERING
vector(abs(vector)>mean(abs(vector(vector~=0)))*outl) = 0;
u = real(vector);
v = imag(vector);

% Adaptive Local Median filtering

%kernel = [-1 -1 -1 -1 -1; -1 2 2 2 -1; -1 2 8 2 -1; -1 2 2 2 -1; -1 -1 -1 -1 -1];
kernel = [-1 -1 -1; -1 8 -1; -1 -1 -1];
tmpv = abs(conv2(v,kernel,'same'));
tmpu = abs(conv2(u,kernel,'same'));
Expand All @@ -591,17 +597,20 @@ function pushbutton_start_Callback(hObject, eventdata, handles)
% Let's throw the outlayers out:
u(u_out) = 0; u(v_out) = 0;
v(v_out) = 0; v(u_out) = 0;
vector = u + sqrt(-1)*v;

vector = u + sqrt(-1)*v;
vector=vector'; %addition
res(:,3) = reshape(real(vector),numrows*numcols,1);
res(:,4) = reshape(imag(vector),numrows*numcols,1);
res(:,4) = reshape(-imag(vector),numrows*numcols,1);%modified

% Filtered results will be stored in '.._flt.txt' file
filt_res = res;


vector=vector'; %addition
vector = fill_holes(vector,numrows,numcols);
vector=vector'; %addition
res(:,3) = reshape(real(vector),numrows*numcols,1);
res(:,4) = reshape(imag(vector),numrows*numcols,1);
res(:,4) = reshape(-imag(vector),numrows*numcols,1);%modified


% scale the pixels and apply the dt
Expand All @@ -627,21 +636,24 @@ function pushbutton_start_Callback(hObject, eventdata, handles)
% Save results as ASCII (text) files:
% Final (filtered, interpolated) results
% fid = fopen([dirname,filesep,filenames(fileind,1:end-4),baseext],'w');

basename = handles.files{fileind}(1:end-4);
baseext = '.vec';

final = fullfile(handles.path,[basename,baseext]);
finaltxt = fullfile(handles.path,[basename,'_final.txt']);
write_openpiv_vec(final,res,xUnits,tUnits,numrows,numcols);
write_openpiv_txt(finaltxt,res,xUnits,tUnits,numrows,numcols);

% Unfiltered, uninterpolated: (comment with % sign if you don't need it)
nofilt = fullfile(handles.path,[basename,'_noflt.txt']);
write_openpiv_vec(nofilt,no_filt_res,xUnits,tUnits,numrows,numcols);
write_openpiv_txt(nofilt,no_filt_res,xUnits,tUnits,numrows,numcols);


% Filtered, but not interpolated:
filtered = fullfile(handles.path,[basename,'_flt.txt']);
write_openpiv_vec(filtered,filt_res,xUnits,tUnits,numrows,numcols);
%filtered2 = fullfile(handles.path,[basename,'_flt2.txt']);
write_openpiv_txt(filtered,filt_res,xUnits,tUnits,numrows,numcols);
%write_openpiv_txt(filtered2,filt_res2,xUnits,tUnits,numrows,numcols);


% Results visualization
Expand Down Expand Up @@ -1207,7 +1219,7 @@ function next_image_Callback(hObject, eventdata, handles)


% --------------------------------------------------------------------
function exit_Callback(hObject, eventdata, handles)
function exit_Callback(~, ~, handles)
% hObject handle to exit (see GCBO)
% eventdata reserved - to be defined in a future version of MATLAB
% handles structure with handles and user data (see GUIDATA)
Expand Down
46 changes: 46 additions & 0 deletions write_openpiv_txt.m
Original file line number Diff line number Diff line change
@@ -0,0 +1,46 @@
function success = write_openpiv_txt(filename,data,xUnits,tUnits,numrows,numcols)
% success = write_openpiv_txt(filename,data,,xUnits,uUnits)
% originally it was just
% fid = fopen(fullfile(handles.path,[handles.files{fileind}(1:end-4),'.txt']),'w');
% fprintf(fid,'%3d %3d %7.4f %7.4f %7.4f\n',res');
% fclose(fid);

% minimal implementation
% filename is a full path
% fullfile(handles.path,[handles.files{fileind}(1:end-4),'.txt'])
fid = fopen(filename,'w');

% data is called res in OpenPIV and it shall be N rows x 5 cols
% in print, it's rotated to 5 rows x N cols
if size(data,2) == 6
data = data';
elseif size(data,1) ~= 6
error('Wrong number of columns');
end

% print the header
% example of the VEC header
% TITLE="E:\2CM_FP500_5%G_68K\C001H001S0015CC\Soapfilmone\Analysis\Run00000
% 1.T000.D000.P000.H001.L.vec" VARIABLES="X mm", "Y mm", "U m/s", "V m/s",
% "CHC", DATASETAUXDATA Application="PIV" DATASETAUXDATA
% SourceImageWidth="1024" DATASETAUXDATA SourceImageHeight="1024"
% DATASETAUXDATA MicrometersPerPixelX="19.530001" DATASETAUXDATA
% MicrometersPerPixelY="19.530001" DATASETAUXDATA LengthUnit="mm"
% DATASETAUXDATA OriginInImageX="0.000000" DATASETAUXDATA
% OriginInImageY="0.000000" DATASETAUXDATA
% MicrosecondsPerDeltaT="2000.000000" DATASETAUXDATA TimeUnit="ms"
% DATASETAUXDATA SecondaryPeakNumber="0" DATASETAUXDATA
% DewarpedImageSource="0" ZONE I=63, J=63, F=POINT

zone = sprintf('ZONE I=%d, J=%d',numrows,numcols);
header = sprintf('VARIABLES= "X %s", "Y %s", "U %s/%s", "V %s/%s", "corr", "CHC", %s\n',...
xUnits,xUnits,xUnits,tUnits,xUnits,tUnits,zone);

% "CHC"

%fprintf(fid,'%s',header);
% print the data
fprintf(fid,'%3d %3d %7.4f %7.4f %7.4f %7.4f\n',data);
fclose(fid);

if fid ~= -1, success = true; else success = false; end
8 changes: 4 additions & 4 deletions write_openpiv_vec.m
Original file line number Diff line number Diff line change
Expand Up @@ -12,9 +12,9 @@

% data is called res in OpenPIV and it shall be N rows x 5 cols
% in print, it's rotated to 5 rows x N cols
if size(data,2) == 5
if size(data,2) == 6
data = data';
elseif size(data,1) ~= 5
elseif size(data,1) ~= 6
error('Wrong number of columns');
end

Expand All @@ -33,14 +33,14 @@
% DewarpedImageSource="0" ZONE I=63, J=63, F=POINT

zone = sprintf('ZONE I=%d, J=%d',numrows,numcols);
header = sprintf('VARIABLES= "X %s", "Y %s", "U %s/%s", "V %s/%s", "CHC", %s\n',...
header = sprintf('VARIABLES= "X %s", "Y %s", "U %s/%s", "V %s/%s", "corr", "CHC", %s\n',...
xUnits,xUnits,xUnits,tUnits,xUnits,tUnits,zone);

% "CHC"

fprintf(fid,'%s',header);
% print the data
fprintf(fid,'%3d %3d %7.4f %7.4f %7.4f\n',data);
fprintf(fid,'%3d %3d %7.4f %7.4f %7.4f %7.4f\n',data);
fclose(fid);

if fid ~= -1, success = true; else success = false; end