Error While running a while loop (index exceeds number of array elements)

I am using a while loop to remove points from a matrix base on their proximity to another group of points in a separate matrix. However, I get an error at (different iteration numbers sometimes, but usually when itr = 9) that says my index exceeds the number of array elements, but that shouldn't be possible since the size of the subsequent array is determined by the number of groups that the splitapply command creates beforehand. I tried forcibly resetting all the indexing variables to 0 at the end of each loop to see if that would work but it did not. I think it has to something to do with the vectors repeating and thus messing up the indexing, but Im not sure. The error i get is below along with my code. Any help would be greatly appreciated.

Index exceeds the number of array elements (27).

Error in (line 307)
    newFN_prime2(i,:) = FN_BNScombined2(i,:) + flow_node_array(GroupList2(i),:);
    
% Import Data txt file 
BioM_3D_250076_1 = importdata('BioM_3D_250076_gen1.txt');
x = BioM_3D_250076_1.data(:, 1);
y = BioM_3D_250076_1.data(:, 2);
z = BioM_3D_250076_1.data(:, 3);

% Create Hemisphere Domain

[x_dom,y_dom,z_dom] = sphere(80);   % Create Sphere
x_dom = x_dom(41:end,:);        % Keep top 41 x points
y_dom = y_dom(41:end,:);        % Keep top 41 y points
z_dom = z_dom(41:end,:);        % Keep top 41 z points
hemisphere_radius = 80;

figure();
Hemi_sf = surf(hemisphere_radius.*x_dom,hemisphere_radius.*y_dom,hemisphere_radius.*z_dom, 'FaceColor','#4DBEEE','EdgeColor', 'none');
alpha 0.2

axis equal
x_ax_lab = xlabel('x axis', 'Color', '#4DBEEE');
y_ax_lab = ylabel('y axis', 'Color', '#4DBEEE');
z_ax_lab = zlabel('z axis', 'Color', '#4DBEEE');

% Plot Outerboundary Circular Plane
x_c = 0;
y_c = 0;
z_c = 0;
radii_plane = 80;
center_plane = [x_c, y_c]; % center point of circular plane
viscircles(center_plane, radii_plane, 'color', '#77AC30'); 

hold on

Center_Root = scatter3(x_c, y_c, z_c, 'filled', 'k');
%%
% Create Set of Random Points to be Plotted inside Hemisphere

n = 2500;
A = zeros(n,3);

% Generate random points  
%radius = 76;
%s = randn(n,3); 
%r = (rand(n,1)).^(1/3);
%r = r./sqrt(sum(s.^2,2));
%r = radius*(s.*r);

% Extract the x, y, and z coordinates from the array.
%x = r(:,1); % Extract x from column #1.
%y = r(:,2); % Extract y from column #2.
%z = abs(r(:,3)); % Extract z from column #3.

A = [x, y, z]; % stores random points generated 
BNS_gen1 = A;
BioM_3Dgraph1 = scatter3(x, y, z, 'filled', 'r');
 %%
 
%BioM_3D_250076_gen1 = table(x, y, z, 'VariableNames', {'x', 'y', 'z'});
%writetable(BioM_3D_250076_gen1, 'BioM_3D_250076_gen1.txt');
 %%
elim_dist1a = nan(numel(x)); % places nans on the diagonal after distance calculation
BNS_BNS_proximity = 1.5;
for i = 1:n  % looping from largest index to avoid calculating the size of the elim_dist2 matrix without pying the price of dynamic growing of a matrix
  for j = 1:(i-1) 
    elim_dist1a(i,j) = sqrt((x(i)-x(j)).^2 + (y(i)-y(j)).^2 + (z(i)-z(j)).^2);
    elim_dist1a(j,i) = elim_dist1a(i,j);
  end
end
% find the points that have its nearest neighbour further away than HS_HS_threshold:
i2keepBNS_BNS = min(elim_dist1a) > BNS_BNS_proximity;
% put those into one pair of arrays
keep_x1 = x(i2keepBNS_BNS);
keep_y1 = y(i2keepBNS_BNS);
keep_z1 = z(i2keepBNS_BNS);
% and the others into another pair of arrays
x_close_neighbors1 = x(~i2keepBNS_BNS);
y_close_neighbors1 = y(~i2keepBNS_BNS);
z_close_neighbors1 = z(~i2keepBNS_BNS); 
B = [keep_x1, keep_y1, keep_z1];
SizeofB = size(B,1);
%%
% Retains one BNS from each of the HS pairs that fail proximity test 

proximity_matrix = triu(elim_dist1a);
proximity_matrix(proximity_matrix == 0) = nan;
[proximity_rowx, proximity_coly] = find(proximity_matrix < BNS_BNS_proximity);
BNS2replace1 = [proximity_rowx, proximity_coly];
%%
% Removes all duplicated elements from indexing vector

if BNS2replace1 ~= 0
BNS2replace_indx = BNS2replace1; % retains BNS2replace matrix for error tracking
BNS_distinct = unique(BNS2replace_indx(:)); % grabs unique values of BNS2replace_indx
BNS2replace_indx(ismember(BNS2replace_indx,BNS_distinct(histc(BNS2replace_indx(:),BNS_distinct)>1))) = Inf; % changes repeated values to infinity
BNS2return_index = min(BNS2replace_indx,[], 2); % selects min index values of each row
BNS2return_index = BNS2return_index(isfinite(BNS2return_index)); % removes Infinite values from return index vector
sizeBNS2return_index = size(BNS2return_index, 1);

for k = 1:sizeBNS2return_index
    BNS2return1(k,:) = A(BNS2return_index(k),:);
end
end
%%
% adds BNS2Return back to B
if BNS2replace1 ~= 0
B = [B;BNS2return1];
keep_x1 = B(:,1); 
keep_y1 = B(:,2);
keep_z1 = B(:,3);
end
%%
% delete the old BNS plot 

delete(BioM_3Dgraph1); 
%%
% Plot Retained Bubble Nucleation Sites (BNS)
 
BioM_3Dgraph1 = scatter3(keep_x1, keep_y1, keep_z1, 'filled', 'r');
SizeofB = size(B,1);
%% 
% Creates array for storing all the flow nodes

flow_node_array = [x_c, y_c, z_c];
FNA_x = flow_node_array(:,1);
FNA_y = flow_node_array(:,2);
FNA_z = flow_node_array(:,3);
SizeofFNA = size(flow_node_array,1);
%%
% FN-BNS Elimination BNS Gen 1 Test (all Flow nodes vs all BNS)

% Circle of influence elimination FN-BNS test & isolates BNS that fail condition
elim_dist2a = nan(numel(keep_x1));
FN_BNS_threshold = 2;
for i = 1:SizeofB % looping from largest index lets us avoid calculating the size of the elim_dist2 matrix without pying the price of dynamic growing of a matrix
  for j = 1:SizeofFNA
     elim_dist2a(j,i) = sqrt((keep_x1(i)-FNA_x(j)).^2 + (keep_y1(i)-FNA_y(j)).^2 +(keep_z1(i)-FNA_z(j)).^2);  
  end
end
% find the points that have its nearest neighbour further away than FN_BNS_threshold:
i2keepFN_BNS1 = min(elim_dist2a)> FN_BNS_threshold;
% put those into one pair of arrays
keep_x2 = keep_x1(i2keepFN_BNS1);
keep_y2 = keep_y1(i2keepFN_BNS1); 
keep_z2 = keep_z1(i2keepFN_BNS1);
% and the others into another pair of arrays
x_close_neighbors2 = keep_x1(~i2keepFN_BNS1);
y_close_neighbors2 = keep_y1(~i2keepFN_BNS1);
z_close_neighbors2 = keep_z1(~i2keepFN_BNS1);
C = [keep_x2, keep_y2, keep_z2]; 
BNS_kept = C; 
%%
% remove BNS that fail retention criteria

delete(BioM_3Dgraph1);
%%
% Replot retained Bubble Nucleation Sites (BNS)

BioM_3Dgraph1 = scatter3(keep_x2, keep_y2, keep_z2, 'filled', 'r');
%%
% Computes distances between each existing BNS and each existing flow sphere 
 
distFNA_BNS1 = pdist2(flow_node_array, BNS_kept); 
%%
% Split HS_Kept into groups that are nearest to vein nodes in Vein Node Array

[~, minRowIdx] = min(distFNA_BNS1,[],1); 
[GroupID, GroupList1] = findgroups(minRowIdx); 
FN_BNSNeighborGroups1 = splitapply(@(x){x},BNS_kept,GroupID(:));
SizeFN_BNSG1 = size(FN_BNSNeighborGroups1,1);
%% 
% Calculate vector addition and combined normalized vectors for each VN-HS group  

for i = 1:SizeFN_BNSG1 
    FN_BNSvectors1{i,:} = FN_BNSNeighborGroups1{i}-flow_node_array(GroupList1(i),:); 
    FN_BNSaddition1(i,:) = sum(FN_BNSvectors1{i},1);
end
%%
% Creates combined normal vectors for new vein nodes
 
FN_BNScombined1 = 1*(FN_BNSaddition1./vecnorm(FN_BNSaddition1,2,2));
SizeFN_BNS_Comb1 = size(FN_BNScombined1, 1);
%%
% Adds new vein node unit vectors to the vein node they are calculated from

for i = 1:SizeFN_BNS_Comb1
    newFN_prime1(i,:) = FN_BNScombined1(i,:) + flow_node_array(GroupList1(i),:);
end
%%
% test for finding repeat combined flow nodes/vectors

if any(ismember(newFN_prime1,flow_node_array)) % if any values in the 2 arrays are equal this code should activate the loop
   FN_FN_threshold = 0.25; % if distance between VNs is greater than this value then point should be retained 
   SizeNewFN_prime = size(newFN_prime1, 1); 
   newFN_primex = newFN_prime1(:,1);
   newFN_primey = newFN_prime1(:,2);
   newFN_primez = newFN_prime1(:,3);
   reconfig_dist = nan(numel(FNA_x)); % preallocation step
for j = 1:SizeofFNA 
  for i = 1:SizeNewFN_prime
     reconfig_dist(j,i) = sqrt((newFN_primex(i)-FNA_x(j)).^2 + (newFN_primey(i)-FNA_y(j)).^2 + (newFN_primez(i)-FNA_z(j)).^2);   
  end  
end
% find the points that are greater than FN_FN_threshold
i2keepFN_FN = find(min(reconfig_dist)> FN_FN_threshold); % keeps the FN points that pass the test 
inot2keepFN_FN = find(min(reconfig_dist)< FN_FN_threshold); % keeps the FN points that fail the test 
% store points greater than FN_FN_threshold
x_newFN_retained = newFN_primex(i2keepFN_FN);
y_newFN_retained = newFN_primey(i2keepFN_FN); 
z_newFN_retained = newFN_primez(i2keepFN_FN); 
% store points less than FN_FN_threshold
x_notkept = newFN_primex(~i2keepFN_FN);
y_notkept = newFN_primey(~i2keepFN_FN);
z_notkept = newFN_primez(i2keepFN_FN);
newFN_retained = [x_newFN_retained, y_newFN_retained, z_newFN_retained]; % store new vein nodes the pass test

FN_VectorGroupNum = inot2keepFN_FN; % FN group # needs to be kept for linking back to proper VNA point

FN_reconfig = FN_BNSvectors1{FN_VectorGroupNum}; % grabs the vector group that makes up repeated combined vector.
FN_BNScombBeta = 1*(FN_reconfig./vecnorm(FN_reconfig,2,2)); % turns vector group into scaled unit vectors
SizeFN_BNScombBeta = size(FN_BNScombBeta,1);

FNconfig_beta = FN_BNScombBeta(:,:) + flow_node_array(GroupList1(FN_VectorGroupNum),:); % adds reconfigured vectors to appropriate vein node chain

% Reconstructs newVN_prime 
newFN_prime1 = [newFN_retained;FNconfig_beta];
SizeFN_BNS_Comb1 = size(newFN_prime1,1);
end % end of if statement 
%%
% Plots new vein node coordinates and their corresponding circles

scatter3(newFN_prime1(:, 1), newFN_prime1(:, 2), newFN_prime1(:, 3), 'filled', 'k');
%%
% appends new vein nodes to array for storing all the vein nodes

flow_node_array = [flow_node_array; newFN_prime1;];
FNA_x = flow_node_array(:, 1);
FNA_y = flow_node_array(:, 2);
FNA_z = flow_node_array(:, 3);
SizeofFNA = size(flow_node_array,1);
newFN_prime_x = newFN_prime1(:, 1);
newFN_prime_y = newFN_prime1(:, 2);
newFN_prime_z = newFN_prime1(:, 3);
SizeNewFN_prime = size(newFN_prime1,1);
%%
% Circle of influence elimination test between newest FN(s) and all BNS & isolates hormone seeds that fail condition
SizeofC = size(C,1); 
elim_dist3a = nan(numel(keep_x2));
FN_BNS_threshold = 2;
for i = 1:SizeofC % looping from largest index lets us avoid calculating the size of the elim_dist2 matrix without pying the price of dynamic growing of a matrix
  for j = 1:SizeNewFN_prime
     elim_dist3a(j,i) = sqrt((keep_x2(i)-newFN_prime_x(j)).^2 + (keep_y2(i)-newFN_prime_y(j)).^2 + (keep_z2(i)-newFN_prime_z(j)).^2);    
  end
end
% find the points that have its nearest neighbour further away than FN_BNS_threshold:
i2keepFN_BNS2 = min(elim_dist3a)> FN_BNS_threshold;
% put those into one pair of arrays
keep_x3 = keep_x2(i2keepFN_BNS2);
keep_y3 = keep_y2(i2keepFN_BNS2);
keep_z3 = keep_z2(i2keepFN_BNS2); 
% and the others into another pair of arrays
x_close_neighbors3 = keep_x2(~i2keepFN_BNS2);
y_close_neighbors3 = keep_y2(~i2keepFN_BNS2);
z_close_neighbors3 = keep_z2(~i2keepFN_BNS2);
D = [keep_x3, keep_y3, keep_z3]; 
BNS_kept = D; 
%%
% Place X over absorbed Hormone Seeds 

BNSx_absorbed = x_close_neighbors3;
BNSy_absorbed = y_close_neighbors3;
BNSz_absorbed = z_close_neighbors3;
BNS_absorbed = [BNSx_absorbed, BNSy_absorbed]; 
BNSabs_graph = scatter3(BNSx_absorbed, BNSy_absorbed, BNSz_absorbed, 'kX');
%%
SizeBNS_keptCount = size(BNS_kept,1);
itr = 0;
while SizeBNS_keptCount >= 0 

% Computes distances between each existing BNS and each existing flow sphere  
distFNA_BNS2 = pdist2(flow_node_array, BNS_kept); 

% Split BNS_Kept into groups that are nearest to Flow nodes in Flow Node Array

[~, minRowIdx] = min(distFNA_BNS2,[],1); 
[GroupID, GroupList2] = findgroups(minRowIdx); 
FN_BNSNeighborGroups2 = splitapply(@(x){x},BNS_kept,GroupID(:));
SizeFN_BNSG2 = size(FN_BNSNeighborGroups2,1);

% Calculate vector addition and combined normalized vectors for each FN-BNS group  

for i = 1:SizeFN_BNSG2 
    FN_BNSvectors2{i,:} = FN_BNSNeighborGroups2{i}-flow_node_array(GroupList2(i),:); 
    FN_BNSaddition2(i,:) = sum(FN_BNSvectors2{i},1);
end

% Creates combined normal vectors for new flow nodes
 
FN_BNScombined2 = 1*(FN_BNSaddition2./vecnorm(FN_BNSaddition2,2,2));
SizeFN_BNS_Comb2 = size(FN_BNScombined2, 1);

% Adds new flow node unit vectors to the flow node they are calculated from

for i = 1:SizeFN_BNS_Comb2
    newFN_prime2(i,:) = FN_BNScombined2(i,:) + flow_node_array(GroupList2(i),:);
end

% test for finding repeat combined flow nodes/vectors

if any(ismember(newFN_prime2,flow_node_array)) % if any values in the 2 arrays are equal this code should activate the loop
   FN_FN_threshold = 0.25; % if distance between VNs is greater than this value then point should be retained 
   SizeNewFN_prime = size(newFN_prime2, 1); 
   newFN_primex = newFN_prime2(:,1);
   newFN_primey = newFN_prime2(:,2);
   newFN_primez = newFN_prime2(:,3);
   reconfig_dist = nan(numel(FNA_x)); % preallocation step
for j = 1:SizeofFNA 
  for i = 1:SizeNewFN_prime
     reconfig_dist(j,i) = sqrt((newFN_primex(i)-FNA_x(j)).^2 + (newFN_primey(i)-FNA_y(j)).^2 + (newFN_primez(i)-FNA_z(j)).^2);   
  end  
end
% find the points that are greater than FN_FN_threshold
i2keepFN_FN = find(min(reconfig_dist)> FN_FN_threshold); % keeps the FN points that pass the test 
inot2keepFN_FN = find(min(reconfig_dist)< FN_FN_threshold); % keeps the FN points that fail the test 
% store points greater than FN_FN_threshold
x_newFN_retained = newFN_primex(i2keepFN_FN);
y_newFN_retained = newFN_primey(i2keepFN_FN); 
z_newFN_retained = newFN_primez(i2keepFN_FN); 
% store points less than FN_FN_threshold
x_notkept = newFN_primex(~i2keepFN_FN);
y_notkept = newFN_primey(~i2keepFN_FN);
z_notkept = newFN_primez(i2keepFN_FN);
newFN_retained = [x_newFN_retained, y_newFN_retained, z_newFN_retained]; % store new vein nodes the pass test

FN_VectorGroupNum = inot2keepFN_FN; % FN group # needs to be kept for linking back to proper FNA point

FN_reconfig = FN_BNSvectors2{FN_VectorGroupNum}; % grabs the vector group that makes up repeated combined vector.
FN_BNScombBeta = 1*(FN_reconfig./vecnorm(FN_reconfig,2,2)); % turns vector group into scaled unit vectors
SizeFN_BNScombBeta = size(FN_BNScombBeta,1);

FNconfig_beta = FN_BNScombBeta(:,:) + flow_node_array(GroupList2(FN_VectorGroupNum),:); % adds reconfigured vectors to appropriate vein node chain

% Reconstructs newVN_prime 
newFN_prime2 = [newFN_retained;FNconfig_beta];
SizeFN_BNS_Comb2 = size(newFN_prime2,1);
end % end of if statement 

% Plots new flow node(s) coordinates 

scatter3(newFN_prime2(:, 1), newFN_prime2(:, 2), newFN_prime2(:, 3), 'filled', 'k');

% appends new flow nodes to array for storing all the flow nodes

flow_node_array = [flow_node_array; newFN_prime2];
FNA_x = flow_node_array(:, 1);
FNA_y = flow_node_array(:, 2);
FNA_z = flow_node_array(:, 3);
newFN_prime_x = newFN_prime2(:, 1);
newFN_prime_y = newFN_prime2(:, 2);
newFN_prime_z = newFN_prime2(:, 3);
SizeNewFN_prime = size(newFN_prime2,1);
SizeofBNS_kept = size(BNS_kept,1);
BNS_keptx = BNS_kept(:,1);
BNS_kepty = BNS_kept(:,2);
BNS_keptz = BNS_kept(:,3);

% Circle of influence elimination test between newest FN(s) and all BNS & isolates hormone seeds that fail condition 
elim_dist4 = nan(numel(BNS_keptx));
FN_BNS_threshold = 2;
for i = 1:SizeofBNS_kept % looping from largest index lets us avoid calculating the size of the elim_dist2 matrix without pying the price of dynamic growing of a matrix
  for j = 1:SizeNewFN_prime
     elim_dist4(j,i) = sqrt((BNS_keptx(i)-newFN_prime_x(j)).^2 + (BNS_kepty(i)-newFN_prime_y(j)).^2 + (BNS_keptz(i)-newFN_prime_z(j)).^2); 
     
  end
end
% find the points that have its nearest neighbour further away than FN_BNS_threshold:
i2keepFN_BNS2 = min(elim_dist4)> FN_BNS_threshold;
% put those into one pair of arrays
keep_x4 = BNS_keptx(i2keepFN_BNS2);
keep_y4 = BNS_kepty(i2keepFN_BNS2);
keep_z4 = BNS_keptz(i2keepFN_BNS2); 
% and the others into another pair of arrays
x_close_neighbors4 = BNS_keptx(~i2keepFN_BNS2);
y_close_neighbors4 = BNS_kepty(~i2keepFN_BNS2);
z_close_neighbors4 = BNS_keptz(~i2keepFN_BNS2);
E = [keep_x4, keep_y4, keep_z4]; 
BNS_kept = E;
SizeBNS_keptCount = size(BNS_kept,1);
SizeNewFN_prime = 0;
SizeFN_BNS_Comb2 = 0;
SizeFN_BNSG2 = 0;
% Place X over absorbed Hormone Seeds 

BNSx_absorbed = x_close_neighbors4;
BNSy_absorbed = y_close_neighbors4;
BNSz_absorbed = z_close_neighbors4;
BNS_absorbed = [BNSx_absorbed, BNSy_absorbed, BNSz_absorbed]; 
BNSabs_graph = scatter3(BNSx_absorbed, BNSy_absorbed, BNSz_absorbed, 'kX');
itr = itr+1;
end    
    
 

Answer 1

[Vance Blake]

There was a missing line of code that needed to be added to allow comparison test that seaches for repeated vectors to run properly. Added SizeofFNA = size(flow_node_array,1); to calculate the size of the matrix flow_node_array the distance calculation loops to run properly and remove any repeated vectors.

% test for finding repeat combined flow nodes/vectors
if any(ismember(newFN_prime2,flow_node_array)) % if any values in the 2 arrays are equal this code should activate the loop
   FN_FN_threshold = 0.25; % if distance between VNs is greater than this value then point should be retained 
   SizeNewFN_prime = size(newFN_prime2, 1); 
   SizeofFNA = size(flow_node_array,1); % <<<<<<<<<<<<<<<<<Added this line here**********************
   newFN_primex = newFN_prime2(:,1);
   newFN_primey = newFN_prime2(:,2);
   newFN_primez = newFN_prime2(:,3);
   reconfig_dist = nan(numel(FNA_x)); % preallocation step
for j = 1:SizeofFNA 
  for i = 1:SizeNewFN_prime
     reconfig_dist(j,i) = sqrt((newFN_primex(i)-FNA_x(j)).^2 + (newFN_primey(i)-FNA_y(j)).^2 + (newFN_primez(i)-FNA_z(j)).^2);   
  end  
end
% find the points that are greater than FN_FN_threshold
i2keepFN_FN = find(min(reconfig_dist)> FN_FN_threshold); % keeps the FN points that pass the test 
inot2keepFN_FN = find(min(reconfig_dist)< FN_FN_threshold); % keeps the FN points that fail the test 
% store points greater than FN_FN_threshold
x_newFN_retained = newFN_primex(i2keepFN_FN);
y_newFN_retained = newFN_primey(i2keepFN_FN); 
z_newFN_retained = newFN_primez(i2keepFN_FN); 
% store points less than FN_FN_threshold
x_notkept = newFN_primex(~i2keepFN_FN);
y_notkept = newFN_primey(~i2keepFN_FN);
z_notkept = newFN_primez(i2keepFN_FN);
newFN_retained = [x_newFN_retained, y_newFN_retained, z_newFN_retained]; % store new vein nodes the pass test
FN_VectorGroupNum = inot2keepFN_FN; % FN group # needs to be kept for linking back to proper FNA point
FN_reconfig = FN_BNSvectors2{FN_VectorGroupNum}; % grabs the vector group that makes up repeated combined vector.
FN_BNScombBeta = 1*(FN_reconfig./vecnorm(FN_reconfig,2,2)); % turns vector group into scaled unit vectors
SizeFN_BNScombBeta = size(FN_BNScombBeta,1);
FNconfig_beta = FN_BNScombBeta(:,:) + flow_node_array(GroupList2(FN_VectorGroupNum),:); % adds reconfigured vectors to appropriate vein node chain
% Reconstructs newVN_prime 
newFN_prime2 = [newFN_retained;FNconfig_beta];
SizeFN_BNS_Comb2 = size(newFN_prime2,1);
end % end of if statement 

Answer 2

[Alan Stevens]

for j = 1:(i-1) 

when i = 1 this will lead to a j index of zero. Matlab only likes indices of 1 and greater.

[Vance Blake]

Im not sure I follow because that part of the code isnt in the while loop so it has no affect on the loop running

[Alan Stevens]

Ignore my comment - I looked at the wrong part of the code!