How to read and plot multiple files

Hi all,

I have multiple output files that i want to proccess at the same time using loops in sequence. Also i want to plot each results 0.02m apart from each other. I have attached the code and my data. Any help will be appreciated.

[filename, pathname] = uigetfile('*.out', 'Select output file to plot');
fullfilename = strcat(pathname, filename);
if filename ~= 0
    header.title = h5readatt(fullfilename, '/', 'Title');
    header.iterations = double(h5readatt(fullfilename,'/', 'Iterations'));
    tmp = h5readatt(fullfilename, '/', 'dx_dy_dz');
    header.dx = tmp(1);
    header.dy = tmp(2);
    header.dz = tmp(3);
    header.dt = h5readatt(fullfilename, '/', 'dt');
    header.nsrc = h5readatt(fullfilename, '/', 'nsrc');
    header.nrx = h5readatt(fullfilename, '/', 'nrx');

    % Time vector for plotting
    time = linspace(0, (header.iterations - 1) * header.dt, header.iterations)';

    % Initialise structure for field arrays
    fields.ex = zeros(header.iterations, header.nrx);
    fields.ey = zeros(header.iterations, header.nrx);
    fields.ez = zeros(header.iterations, header.nrx);
    fields.hx = zeros(header.iterations, header.nrx);
    fields.hy = zeros(header.iterations, header.nrx);
    fields.hz = zeros(header.iterations, header.nrx);
  for n=1:header.nrx
        path = strcat('/rxs/rx', num2str(n));
        tmp = h5readatt(fullfilename, path, 'Position');
        header.rx(n) = tmp(1);
        header.ry(n) = tmp(2);
        header.rz(n) = tmp(3);
        path = strcat(path, '/');
        fields.ex(:,n) = h5read(fullfilename, strcat(path, 'Ex'));
        fields.ey(:,n) = h5read(fullfilename, strcat(path, 'Ey'));
        fields.ez(:,n) = h5read(fullfilename, strcat(path, 'Ez'));
        fields.hx(:,n) = h5read(fullfilename, strcat(path, 'Hx'));
        fields.hy(:,n) = h5read(fullfilename, strcat(path, 'Hy'));
        fields.hz(:,n) = h5read(fullfilename, strcat(path, 'Hz'));
 fh1=figure('Name', strcat('rx', num2str(n)));
 plot(time, fields.ez(:,n), 'r', 'LineWidth', 2), grid on, xlabel('Time [s]'), ylabel('Field strength [V/m]'), title('E_z')
 set(ax,'FontSize', 16, 'xlim', [0 time(end)]);
  end
end

[Simon Chan]

Do you want to plot all results on one figure? And your xlabel mentioned the data is time (not distance) in x-axis, how do you want the data are plotted 0.02m apart?

[Kabit Kishore]

Hi I want to plot the time in y axis and separate each plot by 0.02m on x axis. As every measurements were taken 0.02m interval.

[Simon Chan]

Your data contains information about Field strength [V/m] and time [s] only. Do you need three axis on one plot with a third axis with distance [m]?

[Kabit Kishore]

Hi Simon. Yes and also I want to carry out this process using loops since I have lots of output files

Answer 1

[Simon Chan]

Try the following:

You may now select multiple files at one time.

[multifilename, pathname] = uigetfile('*.out', 'Select output file to plot','MultiSelect', 'on');
allfilename = strcat(pathname, multifilename);
Nz = length(allfilename);
if ~isempty(allfilename)
    for k = 1:Nz
        fullfilename = allfilename{k};
        header.title = h5readatt(fullfilename, '/', 'Title');
        header.iterations = double(h5readatt(fullfilename,'/', 'Iterations'));
        tmp = h5readatt(fullfilename, '/', 'dx_dy_dz');
        header.dx = tmp(1);
        header.dy = tmp(2);
        header.dz = tmp(3);
        header.dt = h5readatt(fullfilename, '/', 'dt');
        header.nsrc = h5readatt(fullfilename, '/', 'nsrc');
        header.nrx = h5readatt(fullfilename, '/', 'nrx');
        % Time vector for plotting
        time = linspace(0, (header.iterations - 1) * header.dt, header.iterations)';
        % Initialise structure for field arrays
        fields.ex = zeros(header.iterations, header.nrx);
        fields.ey = zeros(header.iterations, header.nrx);
        fields.ez = zeros(header.iterations, header.nrx);
        fields.hx = zeros(header.iterations, header.nrx);
        fields.hy = zeros(header.iterations, header.nrx);
        fields.hz = zeros(header.iterations, header.nrx);
        distance = repelem(0.02*(k-1),length(time),1);                  % Added this line
        for n=1:header.nrx
            path = strcat('/rxs/rx', num2str(n));
            tmp = h5readatt(fullfilename, path, 'Position');
            header.rx(n) = tmp(1);
            header.ry(n) = tmp(2);
            header.rz(n) = tmp(3);
            path = strcat(path, '/');
            fields.ex(:,n) = h5read(fullfilename, strcat(path, 'Ex'));
            fields.ey(:,n) = h5read(fullfilename, strcat(path, 'Ey'));
            fields.ez(:,n) = h5read(fullfilename, strcat(path, 'Ez'));
            fields.hx(:,n) = h5read(fullfilename, strcat(path, 'Hx'));
            fields.hy(:,n) = h5read(fullfilename, strcat(path, 'Hy'));
            fields.hz(:,n) = h5read(fullfilename, strcat(path, 'Hz'));
            plot3(distance,fields.ez(:,n),time);                        % Use function plot3
            grid on;
            hold on;
        end
    end
    xlabel('Distance [m]');
    zlabel('Time [s]');
    ylabel('Field strength [V/m]');
    title('E_z')
    ax = gca;
    set(ax,'FontSize', 12, 'zlim', [0 time(end)]);
end

[Kabit Kishore]

Thank you it works perfectly. I just want to add one more step to this. I want to determine the peaks for each of this signal and want to now at what specific time the peak occurs using 'findpeak' method and plot ontop of the signal. Any idea how this can be done.

[Simon Chan]

I don't have the signal processing toolbox and hence cannot use function findpeak. However, I use function max to simulate similar thing as follows:

%
[peakvalue,peakpos] = max(fields.ez(:,n));                      % Calculate max value and its position
plot3(distance,fields.ez(:,n),time);
grid on;
hold on;
plot3(distance,peakvalue,time(peakpos),'r*','MarkerSize',12);   % Add the maximum point
%

[Kabit Kishore]

Thank you very much

[Kabit Kishore]

Sorry for asking again. But is there a way we can determine or locate small peaks that are appearing in the signal.

[Simon Chan]

I need to use function islocalmax this time as follows:

idx.maxpts = islocalmax(fields.ez(:,n),'MaxNumExtrema',2);
plot3(distance,fields.ez(:,n),time);
grid on;
hold on;
plot3(repelem(distance,1,2),fields.ez(idx.maxpts),time(idx.maxpts),'r*','MarkerSize',12);

[Kabit Kishore]

Thank you for helping. I am still struggling to detect smaller peaks as shown in the picture attached.

[Voss]

try:

idx.maxpts = islocalmax(fields.ez(:,n),'MaxNumExtrema',4);

[Kabit Kishore]

Hi Benjamin. I have tried the above command however it is showing me error as below.

Error using plot3

Vectors must be the same length.

Error in untitled3 (line 46)

plot3(repelem(distance,1,2),fields.ez(idx.maxpts),time(idx.maxpts),'r*','MarkerSize',12);

[Voss]

try this:

plot3(repelem(distance,1,nnz(idx.maxpts)),fields.ez(idx.maxpts),time(idx.maxpts),'r*','MarkerSize',12);

[Kabit Kishore]

Thank you Benjamin. It works perfectly.

[Kabit Kishore]

Hi, i want to create a for loop within a loop where i want to use all the time(idx.maxpts) points and create a new plot. I have tried adding the code as highlighed by bold colour however it massess up the other plots in the loop and only plots a singal point.

plot3(distance,fields.ez(:,n),time); % Use function plot3

grid on;

hold on;

idx.maxpts = islocalmax(fields.ez(:,n),'MaxNumExtrema',4);

for ii=(time(idx.maxpts))

b=time(idx.maxpts)

g=b(2,:)

a=b(1,:)

e=((g-a +0.04)/0.04)^2

d=0.4;

H=0.9;

plot3(d,e,H,'*','Marker','hexagram')

end

grid on;

hold on;

plot3(repelem(distance,1,nnz(idx.maxpts)),fields.ez(idx.maxpts),time(idx.maxpts),'r*','MarkerSize',12);

[Simon Chan]

What is your expected output on the figure of your additional code?

Your previous code has x,y and z axis based on distance, field strength and time.

With your additional code, all 3 axes are changed to d (=0.4, a constant), e(some calculation based on time) and H (=0.9, another constant). So obviously the output is going to be messes up if you want to put everything on the same figure.

Better to verify again or may be you can plot the new data on another figure to see the result first.

[Kabit Kishore]

Hi Simon, i am trying to plot the results on a new figure window. However, when i run the program i get only the last values from the loop. I think there is some problem with using hold on command. Also i want to plot distance and e value at height H=0.9.