- Author:
- aram148 <42922407+aram148@users.noreply.github.com>
- Date:
- 2022-07-22 15:47:05+12:00
- Desc:
- Added documentation for VSM model
- Permanent Source URI:
- https://models.physiomeproject.org/workspace/6b0/rawfile/a8a92308e217ac5626809237dd90a31240b22834/peripheral_airway/Peripheral_matlab/createConnMatrices.m
% A function that will compute the necessary connectivity matrices of a
% symmetric tree given its order, and will output all 5 connectivity
% matrices.
% This function will be structured as a recursive function, with the base
% case being the trivial 3-branch case (order 2)
% DATE: 18/01/18
% AUTHOR: Austin Ibarra
% VERSION: 1.1
% == LOG CHANGE == %
% v1.1: - Updated the code because we now have Horsfield constants up to
% order 16
function [Cjplus,Cjminus,Ciplus,Ciminus1,Ciminus2] = createConnMatrices(order)
%% Create the Connectivity Matrices based on order
if order < 2
%% Throw an error if the depth of the tree is less than 2
error('A tree must be at least of order 2!');
elseif order > 16
%% Throw an error if the depth of the tree is greater than 16
error('Currently, a tree cannot be greater than order 16, because the Horsfield constants for higher orders have yet to be coded in.');
elseif order == 2
%% The base case
Cjplus = [1 1 0]';
Cjminus = [0 0 1]';
Ciplus = [0 0 1];
Ciminus1 = [1 0 0];
Ciminus2 = [0 1 0];
else
%% Call the function again recursively
M = 2^(order-1) - 1; % Number of junctions
N = 2^order - 1; % Number of airways
[A,~,~,F,~] = createConnMatrices(order - 1);
B = zeros((N+1)/2,(M-1)/2);
C = zeros(N,1);
C((N-1)/2) = 1;
C((N+1)/2) = 1;
D = zeros((N+1)/2,M);
E = eye((N-1)/2);
G = zeros(1,N);
G((N-1)/2) = 1;
%% Concatenating the appropriate matrices
ApB = cat(1,A,B);
BpA = cat(1,B,A);
Cjplus = cat(2,ApB,C,BpA);
Cjminus = cat(1,D,E);
Ciplus = Cjminus';
FpBt = cat(2,F,B');
BtpF = cat(2,B',F);
Ciminus1 = cat(1,FpBt,G,BtpF);
Ciminus2 = circshift(Ciminus1,1,2);
end
