From ed9eed5a44dbe8f336ee143c1db4c94de40b39a5 Mon Sep 17 00:00:00 2001
From: David Nolte <nolte@purdue.edu>
Date: Fri, 5 Feb 2021 15:11:42 -0500
Subject: [PATCH] Network SIRS updates

---
 drawnet.m     | 187 ++++++++++++++++++++++++++++++++++++++++++++++++++
 heaviside0.m  |  22 ++++++
 heaviside1.m  |  18 +++++
 maskbilevel.m |  13 ++++
 randbin2.m    |  18 +++++
 sublink.m     |  63 +++++++++++++++++
 subnode.m     |  22 ++++++
 7 files changed, 343 insertions(+)
 create mode 100644 drawnet.m
 create mode 100644 heaviside0.m
 create mode 100644 heaviside1.m
 create mode 100644 maskbilevel.m
 create mode 100644 randbin2.m
 create mode 100644 sublink.m
 create mode 100644 subnode.m

diff --git a/drawnet.m b/drawnet.m
new file mode 100644
index 0000000..ab29c85
--- /dev/null
+++ b/drawnet.m
@@ -0,0 +1,187 @@
+% function y = drawnet(node,opt)
+% draws a graph of the node network
+% opt = 0 draws all nodes in a circle (default if no option stated)
+% opt = 1 draws nodes at radius proportional to degree
+% opt = 2 draw circle of nodes with color according to node.value
+% See also:
+%    DynamicDrawnet.m
+%    drawnetC.m
+
+function drawnet(node,opt)
+
+if nargin < 2
+    opt = 0;
+end
+
+[dum,N] = size(node);
+
+% if opt == 2
+%     mx = node(1).value;
+%     for loop = 2:N
+%         if node(loop).value > mx
+%             mx = node(loop).value;
+%         end
+%     end
+%     mn = node(1).value;
+%     for loop = 2:N
+%         if node(loop).value < mn
+%             mn = node(loop).value;
+%         end
+%     end
+% end
+
+
+h = colormap(jet);
+%h = newcolormap('heated');
+fh = gcf;
+clf;
+%figure(fh)
+
+if opt~=2
+    [N,e,avgdegree,maxdegree,mindegree,numclus,meanclus,Lmax,L2,LmaxL2,meandistance] = clusterstats(node);
+    disp(' ')
+    displine('Number of nodes = ',N)
+    disp(strcat('Number of edges = ',num2str(e)))
+    disp(strcat('Mean degree = ',num2str(avgdegree)))
+    displine('Maximum degree = ',maxdegree)
+    disp(strcat('Number of clusters = ',num2str(numclus)))
+    disp(strcat('mean cluster coefficient = ',num2str(meanclus)))
+    disp(strcat('mean distance = ',num2str(meandistance)))
+    disp(' ')
+    disp(strcat('Lmax = ',num2str(Lmax)))
+    disp(strcat('L2 = ',num2str(L2)))
+    disp(strcat('Lmax/L2 = ',num2str(LmaxL2)))
+    disp(' ')
+end
+
+
+if (opt == 0)||(opt == 2)  % Circle of nodes
+    [A,degree,Lap] = adjacency(node);
+    [N,e,avgdegree,maxdegree,mindegree,numclus,meanclus,Lmax,L2,LmaxL2] = clusterstats(node);
+    deltheta = 2*pi/N;
+    
+    % Draw nodes
+    for loop = 1:N
+        theta = deltheta*loop;
+        rad = N;
+        x(loop) = rad*cos(theta);
+        y(loop) = rad*sin(theta);
+
+        if opt == 2
+            colorval = ceil(230*maskbilevel(node(loop).value,0,1,0,1)) + 1;
+        else
+            colorval = ceil(rand*255);
+        end
+
+        crad = pi/1.5;
+        rectangle('Position',[x(loop)-crad,y(loop)-crad,2*crad,2*crad],...
+            'Curvature',[1,1],...
+            'LineWidth',0.1,'LineStyle','-','FaceColor',h(colorval,:))
+    end
+    
+    % Draw links
+    ind = 0;
+     for nloop = 1:N
+        nlink = node(nloop).numlink;
+        deg = degree(nloop);
+        for linkloop = 1:nlink
+            target = node(nloop).link(linkloop);
+            looptarget = target;
+            if deg == maxdegree
+                maxnode(ind+1) = node(nloop).element;
+            end
+            line([x(nloop) x(looptarget)],[y(nloop) y(looptarget)],'Color', 'b');
+        end
+     end
+   
+    axis equal
+
+
+end % end if opt == 0
+
+
+if opt == 1  % Radius proportional to degree
+
+    % Similarity matrix
+    for rowloop = 1:N
+        for coloop = rowloop+1:N
+
+            a = node(rowloop).link;
+            b = node(coloop).link;
+            c = intersect(a,b);
+            [dum,szc] = size(c);
+            D(rowloop,coloop) = szc;
+            D(coloop,rowloop) = szc;
+        end
+    end
+
+
+    [newDis,finallabel] = hierarchcluster(D);
+    [A,degree,Lap] = adjacency(node);
+
+
+    [N,e,avgdegree,maxdegree,mindegree,numclus,meanclus,Lmax,L2,LmaxL2] = clusterstats(node);
+    maxrad = 1.2*2*N/(mindegree+2);
+    %circle(0,0,maxrad);
+
+    % Draw Nodes
+    deltheta = 2*pi/N;
+    for loop = 1:N
+
+        nodelabel(loop) = finallabel(loop);
+        looplabel(nodelabel(loop)) = loop;
+
+        rad = 2*N/(node(nodelabel(loop)).numlink+2);
+        theta = deltheta*loop;
+
+        x(loop) = rad*cos(theta);
+        y(loop) = rad*sin(theta);
+
+        colorval = ceil(rand*255);
+
+        rectangle('Position',[x(loop)-1,y(loop)-1,2,2],...
+            'Curvature',[1,1],...
+            'LineWidth',0.1,'LineStyle','-','FaceColor',h(colorval,:))
+
+        %pause(0.5)
+    end
+
+    % Draw Links
+    for nloop = 1:N
+
+        nlink = node(nodelabel(nloop)).numlink;
+        deg = degree(nodelabel(nloop));
+
+        for linkloop = 1:nlink
+
+            target = node(nodelabel(nloop)).link(linkloop);
+            looptarget = looplabel(target);
+
+            if deg == maxdegree
+                line([x(nloop) x(looptarget)],[y(nloop) y(looptarget)],'Color','r');
+                %nlink
+            else
+                line([x(nloop) x(looptarget)],[y(nloop) y(looptarget)],'Color', 'b');
+
+            end
+
+        end
+
+
+
+    end
+
+    for loop = 1:N
+        theta = deltheta*loop;
+
+        colorval = ceil(rand*255);
+
+        rectangle('Position',[x(loop)-1,y(loop)-1,2,2],...
+            'Curvature',[1,1],...
+            'LineWidth',0.1,'LineStyle','-','FaceColor',h(colorval,:))
+
+    end
+end % if opt == 1
+
+
+
diff --git a/heaviside0.m b/heaviside0.m
new file mode 100644
index 0000000..72e88cb
--- /dev/null
+++ b/heaviside0.m
@@ -0,0 +1,22 @@
+% function y = heaviside0(x)
+% y = 0 at x = 0
+
+function y = heaviside0(x)
+[sy, sx] = size(x);
+
+y = 0;
+if (sx==0)&&(sy==0)
+    y = 0;
+end
+
+for yloop = 1:sy
+    for xloop = 1:sx
+        if x(yloop,xloop)<0
+            y(yloop,xloop) = 0;
+        elseif x(yloop,xloop) == 0
+            y(yloop,xloop) = 0;
+        elseif x(yloop,xloop) > 0
+            y(yloop,xloop) = 1;
+        end
+    end
+end
diff --git a/heaviside1.m b/heaviside1.m
new file mode 100644
index 0000000..b7f39fc
--- /dev/null
+++ b/heaviside1.m
@@ -0,0 +1,18 @@
+% function y = heaviside1(x)
+% y = 1 at x = 0
+
+function y = heaviside1(x)
+[sy, sx] = size(x);
+y = 0;
+for yloop = 1:sy
+    for xloop = 1:sx
+        
+        if x(yloop,xloop)<0
+            y(yloop,xloop) = 0;
+        elseif x(yloop,xloop) == 0
+            y(yloop,xloop) = 1;
+        elseif x(yloop,xloop) > 0
+            y(yloop,xloop) = 1;
+        end
+    end
+end
diff --git a/maskbilevel.m b/maskbilevel.m
new file mode 100644
index 0000000..ce89f03
--- /dev/null
+++ b/maskbilevel.m
@@ -0,0 +1,13 @@
+% maskbilevel.m
+% function [y mask] = maskbilevel(A,lowcut,hicut,vallo,valhi)
+% Assigns values vallo and valhi to outside of cuts
+% mask is binary 1's where inside cutvalues
+
+function [y mask] = maskbilevel(A,lowcut,hicut,vallo,valhi)
+
+% slope = 0.00001*(hicut-lowcut);
+%mtemp = clip(A,lowcut,hicut,slope);
+mtemp = heaviside0(A-lowcut).*heaviside1(hicut-A);
+y = A.*mtemp + valhi*heaviside1(A-hicut).*(1-mtemp) +  vallo*heaviside0(lowcut-A).*(1-mtemp);
+mask = mtemp;
+
diff --git a/randbin2.m b/randbin2.m
new file mode 100644
index 0000000..58b95bc
--- /dev/null
+++ b/randbin2.m
@@ -0,0 +1,18 @@
+% function y = randbin2(M,N,thresh)
+% random binary matrix of length M by N
+% thresh between 0 and 1
+% p(1) = 1-thresh
+% p(0) = thresh
+
+function y = randbin(M,N,thresh)
+
+temp = rand(M,N);
+y = zeros(M,N);
+for mloop = 1:M
+    for nloop = 1:N
+        if temp(mloop,nloop) > thresh
+            y(mloop,nloop) = 1;
+        end
+    end
+end
+
diff --git a/sublink.m b/sublink.m
new file mode 100644
index 0000000..6048b4f
--- /dev/null
+++ b/sublink.m
@@ -0,0 +1,63 @@
+%function node = sublink(node1,node2,node)
+
+function newnode = sublink(node1,node2,node)
+
+templink1 = node(node1).link;
+templink2 = node(node2).link;
+
+tempnode = node;
+clear node
+
+% templink1
+% node2
+% templink2
+% node1
+w1 = where(templink1,node2);
+w2 = where(templink2,node1);
+
+if (w1(1) == -1)|(w2(1) == -1)
+    disp('Error in sublink: node missing')
+    newnode = tempnode;
+    return
+end
+
+nlink1 = tempnode(node1).numlink;
+if nlink1 == 1                          % only one linke in set 
+    tempnode(node1).numlink = 0;
+    tempnode(node1).link = [];
+elseif w1 == nlink1             % target at end of set
+    newnlink = nlink1 - 1;
+    temp1 = templink1(1:newnlink);
+    tempnode(node1).numlink = newnlink;
+    tempnode(node1).link = temp1;
+    tempnode(node1).numlink = newnlink;
+else
+    newnlink = nlink1 - 1;
+    tempnode(node1).link(w1) = tempnode(node1).link(nlink1);
+    temp1 =  tempnode(node1).link(1:newnlink);
+    tempnode(node1).link = temp1;
+    tempnode(node1).numlink = newnlink;
+end
+
+nlink2 = tempnode(node2).numlink;
+if nlink2 == 1                          % only one linke in set 
+    tempnode(node2).numlink = 0;
+    tempnode(node2).link = [];
+elseif w2 == nlink2             % target at end of set
+    newnlink = nlink2 - 1;
+    temp2 = templink2(1:newnlink);
+    tempnode(node2).numlink = newnlink;
+    tempnode(node2).link = temp2;
+    tempnode(node2).numlink = newnlink;
+else
+    newnlink = nlink2 - 1;
+    tempnode(node2).link(w2) = tempnode(node2).link(nlink2);
+    temp2 = tempnode(node2).link(1:newnlink);
+    tempnode(node2).link = temp2;
+    tempnode(node2).numlink = newnlink;
+end
+
+
+newnode = tempnode;
+
+
diff --git a/subnode.m b/subnode.m
new file mode 100644
index 0000000..f67fddb
--- /dev/null
+++ b/subnode.m
@@ -0,0 +1,22 @@
+%function newnode = subnode(nodenum,node)
+
+function newnode = subnode(nodenum,node)
+
+[dum sz] = size(node);
+
+tempnode = node;
+
+for nodeloop = 1:sz
+    nlink = node(nodeloop).numlink;
+    for linkloop = 1:nlink
+        linknum = node(nodeloop).link(linkloop);
+        if (linknum == nodenum)
+            gonode.label(1) = nodenum;
+            gonode.label(2) = node(nodeloop).element; 
+            tempnode = sublink(gonode.label(1),gonode.label(2),tempnode);
+        end
+    end
+end
+
+newnode = tempnode;
+