how to create a mean velocity matrix from multiple .dat files and plot the mean profile?

I have a folder of x no of *.dat files. So far I've imported 1 file and reshaped the columns to get matrices for each column. e.g for 1 dat file:

filename='B00049.dat';

delimiterIn= ' ';

headerlinesIn = 3;

A= importdata(filename,delimiterIn,headerlinesIn);

%% Matrices%%%%

c = 214;%columns i

r = 134;%rows j

dt= 0.0005;

x=reshape(A.data(:,1),r,c);

y=reshape(A.data(:,2),r,c);

u=reshape(A.data(:,3),r,c);

v=reshape(A.data(:,4),r,c);

I now want to import multiple files and creat a mean velocity matrix fro u and v then plot a profile. How can this be done ? Thanks

[jonas]

Can you provide two of your .dat files?

[Ernest Adisi]

thanks for the response, i'm trying to create a path to the folder containing all the files then perform a loop in order to make the 4 matrices for each file, then finally create a mean velocity matrix from the average of the last two columns in the 2 dat files. Hope that makes sense

Answer 1

[jonas]

files=dir('*.dat') %If files are in current dir, otherwise enter entire path

%%Preallocate some cells
x=cell(1,numel(files))
y=cell(1,numel(files))
u=cell(1,numel(files))
v=cell(1,numel(files))

%%Loop over all files
for i=1:numel(files)
  filename=files(i).name
  delimiterIn= ' ';
  headerlinesIn = 3;
  A= importdata(filename,delimiterIn,headerlinesIn);
  c = 214;%columns i
  r = 134;%rows j
  dt= 0.0005;

%%Save results
  x{i}=reshape(A.data(:,1),r,c);
  y{i}=reshape(A.data(:,2),r,c);
  u{i}=reshape(A.data(:,3),r,c);
  v{i}=reshape(A.data(:,4),r,c);
end

When you have all the data stored in cells, you can concatenate along the third dimensinon and take the average. I assume all files share the same size?

V = cat(3, v{:});
mean(V,3)

[Ernest Adisi]

Yes they are. This works brilliantly. Thanks a lot, can I just ask a further question.

1. for V = cat(3, v{:}); has it concatenated the V matrices so like they are stacked like a 'stack of waffles' for example, then calculates a mean matrix ? and the 3 what does it do please? don't understand the matlab explanation of it

[Ernest Adisi]

basically can you explain how the V=... code works please

[Stephen23]

The command cat(3,...) concatenates all of its inputs along the third dimension. It has multiple inputs supplied as a comma-separated list:

https://www.mathworks.com/help/matlab/matlab_prog/comma-separated-lists.html

https://www.mathworks.com/matlabcentral/answers/320713-how-to-operate-on-comma-separated-lists

[Ernest Adisi]

alright thanks I understand that now. lastly, i'm just wondering what the 'v{:}' in the command V = cat(3, v{:}); is telling matlab

[Ernest Adisi]

and also, from the mean matrix can how do i plot a mean profile. i don't if its a contour plot or something else? so it display the different velocity profiles in the regions

[Stephen23]

" i'm just wondering what the 'v{:}' in the command V = cat(3, v{:}); is telling matlab"

See the links in my comment.

[Ernest Adisi]

oh great. Thanks for the help

[jonas]

It's difficult for me to plot because I am not really familiar with your data. Nevertheless, I made an attempt. Is that what you want?

[Xq Yq]=meshgrid(0:6:max(x{1}(:)),0:6:max(y{1}(:)));
Vq=griddata(x{1}(:),y{1}(:),V_mean(:),Xq,Yq);
contourf(Xq,Yq,Vq)

[Ernest Adisi]

Yes, its that, thank you v.much Jonah. i'm just trying to understand the code but that seems similar to what I need to achieve. Thank you

[jonas]

You're welcome. Attempting to understand the code is exactly what you should be doing, so continue with that. The MATLAB documentation is extremely helpful!

The plotting is fairly simple, but looks complex. Just read about meshgrid (creates a grid) and griddata (interpolates scattered data to the grid).

[Ernest Adisi]

%% plot mean velocities %%%% V-mean%% [Xq Yq]=meshgrid(0:2:max(x{1}(:)),0:2:max(y{1}(:))); Vq=griddata(x{1}(:),y{1}(:),V_m(:),Xq,Yq);

f1=figure; contourf(Xq,Yq,Vq); xlabel({'x(m)'}) ylabel({'y(m)'}) title('V-mean velcity profile') hcb= colorbar title(hcb,'Velocity[m/s]')

%%%% U-mean%% [Xq Yq]=meshgrid(0:2:max(x{1}(:)),0:2:max(y{1}(:))); Uq=griddata(x{1}(:),y{1}(:),U_m(:),Xq,Yq);

f2=figure; contourf(Xq,Yq,Uq) xlabel({'x(m)'}) ylabel({'y(m)'}) title('U-mean velcity profile') hcb= colorbar title(hcb,'Velocity[m/s]')

this seems to have a problem with the plot, do you have an idea why please?

[jonas]

Please make sure your code is formated like this:

 [Xq Yq]=meshgrid(0:2:max(x{1}(:)),0:2:max(y{1}(:))); 
 Vq=griddata(x{1}(:),y{1}(:),V_m(:),Xq,Yq);

 f1=figure; 
contourf(Xq,Yq,Vq); 
xlabel({'x(m)'}) ylabel({'y(m)'}) 
title('V-mean velcity profile') 
hcb= colorbar  
title(hcb,'Velocity[m/s]')

%%%%U-mean%%
[Xq Yq]=meshgrid(0:2:max(x{1}(:)),0:2:max(y{1}(:))); 
Uq=griddata(x{1}(:),y{1}(:),U_m(:),Xq,Yq);

f2=figure; 
contourf(Xq,Yq,Uq) 
xlabel({'x(m)'}) 
ylabel({'y(m)'}) 
title('U-mean velcity profile') 
hcb= colorbar 
title(hcb,'Velocity[m/s]')

Seems fine to me, provided that the grid should start at x=0 and y=0. Otherwise, change to:

 [Xq Yq]=meshgrid(min(x{1}(:)):2:max(x{1}(:)),min(x{1}(:)):2:max(y{1}(:)));

Other than that, the code looks fine at first glance. What is the error message?

[Ernest Adisi]

its sorted out now, thanks again

[Ernest Adisi]

Hey, can i ask another question please. I trying to change the orientation of the profile,like 90 degree rotation, could you show me how its done please

[jonas]

What do you mean? Do you want to just change the view? Then change the first value in the view pair. For example:

set(gca,'view',[90 90]) %

...or do you want to rotate the actual Z-matrix 90 degrees?

B = rot90(A)

this will however crash your script unless you have a square matrix. You probably just want to change the view.

[Ernest Adisi]

like as you see in the plot above, the flow direction is top to bottom , but i want it left to right so rotation in z axis. while keeping the x and y the same.

[Ernest Adisi]

i think its the viewpoint functioion but not sure how to go about it

[jonas]

Well I gave you two options! Try changing the view and see if it looks OK

[Ernest Adisi]

for example, this rotated anticloclwise 90 degrees but everything else kept the same

[jonas]

Yes, as I said above.

set(gca,'view',[90 90])

[Ernest Adisi]

set(gca,'view',[90 90] this works but the axis x and y have change places

[Ernest Adisi]

sorry, keep seeing your updated message late, i tried that one you gave me

[jonas]

I dont understand. Do you want the flow direction to be ~0.6-0.7 m in your latest image?

[Ernest Adisi]

you're right. It worked. Thank you so much again for all the help. Appreciate it

[jonas]

No problemo!