MATLAB Script (PIV model)
%% PIV_data_process
% Author: Arif Hossain
clc; close all; clear all;
%% Experimental data
load ('PIV_data.mat')
x_raw = data(:,1);y_raw = data(:,2);Vx_raw = data(:,3);Vy_raw = data(:,4);
x = reshape(x_raw,128,66); ; y = reshape(y_raw,128,66); ; Vx = reshape(Vx_raw,128,66); Vy = reshape(Vy_raw,128,66); V_mag = sqrt(Vx.^2+Vy.^2);
% Contour plot of Velocity magnitude
figure()
contourf(x,y,V_mag,50,'linestyle','none'); hold on
quiver(x(1:4:end),y(1:4:end),Vx(1:4:end),Vy(1:4:end),2,'color',[0 0 0])
set(gca,'fontsize',14);xlabel('x [mm]');ylabel('y [mm]');axis image ; h = colorbar; title(h, 'Velocity Magnitude [m/s]','fontsize',12) ; set(gca,'fontsize',12)
figure() % Contour plot of Vx
contourf(x,y,V_mag,50,'linestyle','none') set(gca,'fontsize',14);xlabel('x [mm]');ylabel('y [mm]');axis image h = colorbar; title(h, 'Vx [m/s]','fontsize',12) set(gca,'fontsize',12)
figure() % Contour plot of Vy
contourf(x,y,Vy,50,'linestyle','none') set(gca,'fontsize',14);xlabel('x [mm]');ylabel('y [mm]');axis image
h = colorbar; title(h, 'Vy [m/s]','fontsize',12) set(gca,'fontsize',12)
%% Image read
A = imread('Image1.bmp'); B = imread('Image2.bmp');
Am = ((A));%(1:707,1:1381,1:3)));
Bm = ((B));%(1:707,1:1381,1:3)));
figure()
subplot (2,1,1);set(gca,'fontsize',12)
image(Am);xlabel('X-pixel');ylabel('Y-pixel');title('Image 1');
xlim([0 1381]);ylim([0 707]);
subplot(2,1,2);set(gca,'fontsize',12)
image(Bm);xlabel('X-pixel');ylabel('Y-pixel');title('Image 2')
xlim([0 1381]);ylim([0 707]);
%% Cross-correlation calculation (32x32,50% overlap)
[yL,xL] = size(Am(1:707,1:1381));
IA = 32; % Interogation area
overlap = 0.5; % 50% overlap
L = 0.3048; % Length
W = 0.1542; % Width
dt= 3.9e-6; % Time
L_pix = L/xL; % Length per pixal
W_pix = W/yL; % Width per pixal
count_y = 1;
for IR_y = 17:8:672 %IA/2+1:IA/4:672
count_x = 1;
for IR_x = 17:8:1376 %IA/2+1:IA/4:1376
IA_1(:,:) = Am(IR_y:IR_y+31,IR_x:IR_x+31);
for IR_yy = 1:IA
for IR_xx = 1:IA
IA_2(:,:) = Bm(IR_y+IR_yy-17:IR_y+IR_yy+14,IR_x+IR_xx-17:IR_x+IR_xx+14);
IA_1_mean = mean(mean(IA_1));
IA_2_mean = mean(mean(IA_2));
IA_1_std = sqrt(sum(sum((IA_1-IA_1_mean).^2)));
IA_2_std = sqrt(sum(sum((IA_2-IA_2_mean).^2)));
R(IR_yy,IR_xx) = sum(sum((IA_1-IA_1_mean).*(IA_2-IA_2_mean)))/(IA_1_std*IA_2_std);
[MaxR,IMax] = max(max(R));
[r,c] = find(R == MaxR);
dx(count_y,count_x) = (IR_x+15)*W_pix;
dy(count_y,count_x) = (yL-(IR_y+15))*W_pix-0.07799;
Vxx(count_y,count_x) = ((c(1)-17)*W_pix)/(3.9*10^(-6));
Vyy(count_y,count_x) = ((r(1)-17)*L_pix)/(3.9*10^(-6));
V_magnitude (count_y,count_x) = sqrt(Vxx(count_y,count_x).^2+Vyy(count_y,count_x).^2);
end
end
count_x = count_x +1;
end
count_y = count_y +1;
end
%% Velocity correction
for i = 1:82
for j = 1:170
if Vxx(i,j)>300
Vxx(i,j) = 0;
elseif Vyy(i,j)> 80
Vyy(i,j) = 0;
elseif V_magnitude(i,j)> 300
V_magnitude(i,j)= 0;
end
end
end
xx = dx*1000; yy = dy*1000;
%% Velocity countor plot
figure()
contourf(xx,yy,V_magnitude,100,'linestyle','none') hold on
quiver(xx,yy,Vxx,Vyy,7,'color',[0 0 0]) set(gca,'fontsize',14);xlabel('x [mm]');ylabel('y [mm]');axis image
h = colorbar; title(h, 'Velocity Magnitude [m/s]','fontsize',12) set(gca,'CLim',[0 300],'fontsize',12); xlim([10 300]); ylim([-70 70])
figure()
contourf(xx,yy,Vxx,100,'linestyle','none'); set(gca,'fontsize',14);xlabel('x [mm]');ylabel('y [mm]');axis image
h = colorbar; title(h, 'Velocity Magnitude [m/s]','fontsize',12); set(gca,'CLim',[0 300],'fontsize',12); xlim([10 300]); ylim([-70 70])
figure()
contourf(xx,yy,Vyy,100,'linestyle','none');set(gca,'fontsize',14);xlabel('x [mm]');ylabel('y [mm]');axis image
h = colorbar; title(h, 'Velocity Magnitude [m/s]','fontsize',12); set(gca,'CLim',[-60 60],'fontsize',12); xlim([10 300]); ylim([-70 70])