% Index-finding routine for a finite difference model of electrical
% impedance tomography (EIT) in a square. We assume that the square is
% divided into NxN pixels with N an odd integer.
%
% Outside edges of boundary pixels, understood as electrodes used in EIT
% measurements are numbered counterclockwise, starting from the middle
% pixel on the right ("Eastern") side of the square. Note that corner
% pixels have only one electrodes.
%
% Given a row index and column index in the NxN grid, this routine finds
% b_ind, the index of the unique electrode attached to the input pixel. 
% If the input pixel is not on the boundary, b_ind=[].
%
% See also: bind2gind.m
%
% Arguments:
% N     grid size parameter: NxN pixels
% row   row index in the NxN grid
% col   column index in the NxN grid
%
% Returns:
% b_ind  number of electrode at the input pixel, range 1,2,3,...,4*(N-1)
%
% Samuli Siltanen March 2025

function b_ind = gind2bind(N,row,col)

% Check that N is odd
if mod(N,2)==0
    error('N should be an odd integer')
end
% Check that row has reasonable value
if (row<1)|(row>N)
    error('Row index should be in the range 1,2,3,...,N')
end
% Check that col has reasonable value
if (col<1)|(col>N)
    error('Column index should be in the range 1,2,3,...,N')
end

% Define halfN so that N = 2*halfN+1
halfN = floor(N/2);

% Initialize result
b_ind = [];

% Left (Western) boundary
if (col==1)&(row>1)
    b_ind = (halfN+1)+(N-1)+(row-1);
end

% Right (Eastern) boundary
if (col==N)&(row<N)
    if row>halfN+1
        b_ind = (halfN+1)+3*(N-1)+(N-row);
    else
        b_ind = halfN+2-row;
    end
end

% Top (Northern) boundary
if (row==1)&(col<N)
    b_ind = (halfN+1)+N-col;
end

% Bottom (Southern) boundary
if (row==N)&(col>1)
    b_ind = (halfN+1)+2*(N-1)+(col-1);
end

