Extracting data from meshgrid output

Question

I am trying to extract data from meshgrid output where zCs is 1000 x 1000 matrix of z values which is produced from meshgrid of s using bin2mat. s is a struct which has x, y, z data. res is 1000.

Based on the threshold I am trying to extract logical which is then used to extract data from s whose z values are between the thresholds.

This code is working but very slow as shown on the screenshot.

xi = linspace(min(s.X),max(s.X),res);
yi = linspace(min(s.Y),max(s.Y),res);  

[~,~, iBins] = histcounts(s.X,xi);
[~,~, jBins] = histcounts(s.Y,yi);

n_pts = length(s.X);

in1 = zCs > 0.0030 & zCs < 0.1500;
Ind = find(in1(:));

[nrows ncols] = size(in1);

I = rem(Ind-1,nrows)+1;
J = (Ind-I)/nrows + 1;

% 
% 
% [I,J] = ind2sub(size(in1), Ind);
N_true = length(I);


hwb = waitbar(0,'Please wait...');

in = false(n_pts,1);
  for k = 1:N_true
      in3 =  iBins == J(k) & jBins == I(k);
      in = in | in3;

      waitbar(k / N_true);

  end
  close(hwb);
  sum(in);




xk = s.X(in);
yk = s.Y(in);
zk = s.Z(in);

Answer 1

Welcome to the world of vectorization!

Your loop checks all elements in iBins and jBins against each of the elements in J and I. You do in = in | in3 to fill in with all the True values.

Now, this is the prime example of something that can be vectorized. Vectorization means that instead of working with one by one element in a vector, you work with the entire vector at the same time. This is much faster in MATLAB, as it's what MATLAB is designed to do.

A very short example (let's disregard jBins and I for a minute):

iBins = [4  2  3  8  6  4  2]
J     = [1  4  2  6  7]

Now, we can loop through J, and do in = in | in3 to get the following (this is the resulting in after each iteration:

k = 1
in3 = 0  0  0  0  0
in  = 0  0  0  0  0
k = 2
in3 = 0  1  0  0  0
in  = 0  1  0  0  0
k = 3
in3 = 0  0  1  0  0
in  = 0  1  1  0  0
k = 4
in3 = 0  0  0  1  0
in  = 0  1  1  1  0
k = 5
in3 = 0  0  0  0  0
in  = 0  1  1  1  0

The last in = 0 1 1 1 0 is the result after the loop.

Instead of looping like this, we can check the vectors against each other, and create a matrix:

    4  2  3  8  6  4  2
    -  -  -  -  -  -  -
1 | 0  0  0  0  0  0  0
4 | 1  0  0  0  0  1  0
2 | 0  1  0  0  0  0  1
6 | 0  0  0  0  1  0  0
7 | 0  0  0  0  0  0  0

If we check each row for true values, we'll see that we have true values in the second, third and fourth row, the same as above. How to do this?

logical_matrix = bsxfun(@eq, iBins, J.');
in = any(logical_matrix, 2);

Note that I transposed J, so that the two arguments to bsxfun are orthogonal.

---

Let's bring back jBins and I:

You must have true values in the same positions in the two matrices, so:

logical_matrix = bsxfun(@eq, iBins, J.') & bsxfun(@eq, jBins, I.');
in = any(logical_matrix, 2);

_{logical_matrix isn't a good name, I just couldn't think of a better word at the moment. You can skip it entirely, and do any(bsxfun..., 2).}

If you have MATLAB R2016b or newer you can simply do:

in = any((iBins == J.') & (jBins == I.'), 2);

instead of bsxfun.

---

The rest of the code:

I recommend having parentheses in the below line, to make it sure what you compare, and make sure MATLAB works in the correct order.

in1 = (zCs > 0.0030) & (zCs < 0.1500);

It behaves the right way in the line above, but it's not always clear which order MATLAB will execute such commands (do == have precedence over < or not?) An example where the order is the other way around:

1 == 1 + 2 == 2
ans = 0
(1 == 1) + (2 == 2)
ans = 2

---

It might make sense when you know what the entire code does, but I had to read the code many times to understand what in and in3 were. The same goes with zCs.

---

You never use ncols, so you can do: nrows = size(in1, 1) or [nrows, ~] = size(in1). This is a very fast operation, so speed doesn't really matter, it's up to you to decide what you think is easiest to read and write. Anyway, MATLAB recommends that you separate output arguments with comma: [nrows, ncols] = size(in1) (I'm quite sure the editor gives you a recommendation (orange line in the scroll bar to the right).

---

waitbar is a time consuming operation, but you might want it there when debugging.