#!/usr/bin/perl

# s-core+ for undirected weighted networks
#
# This is a version of the s-core+ that does not keep a sorted list of the node strengths, but finds the smallest strength by sorting for each s-core. It is slow on large networks with many distinct cores. It keeps a sorted weight list after the innermost s-core is found.
# The version stores A-B and B-A links for simplicity and is memory intensive for large networks.
#
# Inputs:
# 1. Network file path/name
# 2. s-core+ output file path/name
#
# Outputs:
# 1. s-core+ file defined in input 2.
#
use strict;

my $netFile = shift @ARGV;
my $outFile = shift @ARGV;

#-------------------------
# Switches and constants
#
my $ifKCore = "no"; # if "yes", all edges will be set to 1
my $networkSizeForFallback = 6000; # Number of nodes in network needed to make fallback copies needed for s-core+. This should be close to the innermost s-core size as the algorithm runs a lot slower with this.
my $accurateDigits = 9; # number of digits in comparisons
my $ifVerbose = "no"; # if "yes", s-core index and number of remaining nodes will be printed each iteration

my $nodeAColumn = 0; # 0 1 2 etc.
my $nodeBColumn = 1; # 0 1 2 etc.
my $weightColumn = 2;
#-------------------------

my %network;
my %strength;
my %core;
my %thresh;
my %coreSize;
#-------------------------
# s-core plus variables

my %fallbackNetwork;
my %fallbackStrength;

my @sortedWeights = ();

my $threshold;
my $scorePlus = 0; # becomes 1 when s-core+ starts after the pure s-core is done

#-------------------------

# Loading network
&loadNetwork(\%network, \%strength, $netFile, $nodeAColumn, $nodeBColumn, $weightColumn);

#------------
# s-core plus
#------------
my $coreNumber = 0; # s-core is defined so that all nodes are contained in the 0-core, while the 1-core is determined by the weakest node
my $prevCoreSize = scalar(keys(%strength)) + 1;
my $loopCheck = 0;

#----------------------------
# s-core plus main loop start
while () {

	#------------------
	# s-core loop start
	while (scalar(keys(%strength)) > 0) {
		$coreNumber++ if $prevCoreSize != scalar(keys(%strength));
		$prevCoreSize = scalar(keys(%strength));
		$coreSize{$coreNumber} = scalar(keys(%strength));
		if (scalar(keys(%strength)) <= $networkSizeForFallback) {
			&fallback(\%fallbackNetwork, \%fallbackStrength, \%network, \%strength); # fallbacknetwork = network
		}

		# Finding the smallest node strength(s) and storing them in a hash to be used in next step.
		my %checkNodes;
		my $iterator = 0;

		foreach my $node (sort {$strength{$a} <=> $strength{$b}} keys %strength) {

			if ($iterator == 0) {
				$threshold = $strength{$node};
				$checkNodes{$node} = 1;
				$iterator++;
				next;
			}

			if (($iterator != 0) and (&acc($strength{$node}) == &acc($threshold))) {
				$checkNodes{$node} = 1;
				$iterator++;
			} else {
			last;		
			}
		}

		$loopCheck = 0;
		my $prevNumberOfCheckNodes = scalar(keys(%checkNodes)) + 1;
		# iterative removal of nodes with lower node strength than $threshold. These nodes are stored in %checkNodes
		while (scalar(keys(%checkNodes)) > 0) {
			$loopCheck++ if $prevNumberOfCheckNodes == scalar(keys(%checkNodes));
			$prevNumberOfCheckNodes = scalar(keys(%checkNodes));

			# Checking if the loop is looping. If so, it will be forced to end
			if ($loopCheck > 10) {
				print STDERR "\nBEWARE. Nodes remaining: There may be issues with the s-core (plus) results due to accuracy. Printing current results to $outFile. Ending script.\n";
				last;
			}

			foreach my $checkNode (keys %checkNodes) {
				if (not exists $strength{$checkNode}) {
					delete $checkNodes{$checkNode};
					next;
				}
				foreach my $node (keys %{ $network{$checkNode} }) {
					# Here: $checkNode deleted, edges connected to $checkNode deleted, $node stored in %checkNodes if below $threshold
					&updateHashes(\%network, \%strength, \%core, \%thresh, \%checkNodes, $coreNumber, $node, $checkNode, $threshold);
				}
			}
		}
		last if $loopCheck > 10;

		if ($ifVerbose eq "yes") {
		    print STDERR "\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\bindex: " . $coreNumber . "   nodes remaining: " . scalar(keys(%strength)) . "      ";
	    }
	}
	# s-core loop end
	#----------------

	last if $loopCheck > 10;

	if ((scalar(keys(%fallbackStrength)) == 0) and ($scorePlus == 0)) {
		print STDERR "\ns-core finished at s-core index n = $coreNumber, collapsing for a core size of N = ",$prevCoreSize,". \nBEWARE. s-core plus did not run. If you want s-core plus to run, increase \$networkSizeForFallback so that it is larger than ",$prevCoreSize,". The s-core results have been printed to $outFile\n";
		last;
	}

	&fallback(\%network, \%strength, \%fallbackNetwork, \%fallbackStrength); # network = fallbacknetwork

	if ($scorePlus == 0) {
		print STDERR "\ns-core finished at s-core index n = $coreNumber, collapsing for a core size of N = ",$prevCoreSize,". Continuing with s-core plus:\n";
		&sortedLinkWeights(\@sortedWeights, \%network);		
		$scorePlus = 1;
	}

	my $thresholdIsSame = 1;

	while ($thresholdIsSame) {
		my @weakestLink = shift(@sortedWeights);

		# Checking if the nodes connecting the weakest link still exists, if so hashes are updated
		if ((exists $strength{$weakestLink[0][0]}) and (exists $strength{$weakestLink[0][1]})) {
			&updateScorePlusHashes(\%network, \%strength, $weakestLink[0][0], $weakestLink[0][1], $weakestLink[0][2]);

			# Checking whether the link deletion causes a new (smaller) $threshold
			if (&min2($strength{$weakestLink[0][0]}, $strength{$weakestLink[0][1]}) <= $threshold) {
				$thresholdIsSame = 0;
				$threshold = &min2($strength{$weakestLink[0][0]}, $strength{$weakestLink[0][1]});

				# If the next links in @sortedWeights is equal to weakestLink, or nonexisting in %network, they must be removed, and $threshold might be lowered further
				$loopCheck = 0;
				my $prevNumberOfSortedWeights = $#sortedWeights + 1;

				while ((&acc($sortedWeights[0][2]) == &acc($weakestLink[0][2])) or ((not exists $strength{$sortedWeights[0][0]}) or (not exists $strength{$sortedWeights[0][1]}))) {
					$loopCheck++ if $prevNumberOfSortedWeights == $#sortedWeights;
					$prevNumberOfSortedWeights = $#sortedWeights;

					# Checking if the loop is looping. If so, it will be forced to end
					if ($loopCheck > 10) {
						print STDERR "\nBEWARE. Links remaining: There may be issues with the s-core (plus) results due to accuracy. Printing current results to $outFile. Ending script.\n";
						@sortedWeights = ();
						last;
					}
					if ((not exists $strength{$sortedWeights[0][0]}) or (not exists $strength{$sortedWeights[0][1]})) {
						shift(@sortedWeights);
						last if $#sortedWeights < 0;
						next;
					}
					&updateScorePlusHashes(\%network, \%strength, $sortedWeights[0][0], $sortedWeights[0][1], $sortedWeights[0][2]);
					$threshold = &min2($threshold, &min2($strength{$sortedWeights[0][0]}, $strength{$sortedWeights[0][1]}));
					shift(@sortedWeights);
				}
			}
		}
		last if $#sortedWeights < 0;
	}
	last if $#sortedWeights < 0;
}
# s-core plus main loop end
#--------------------------

print STDERR "\n";

# Printing to $outFile
open (OUT, '>' . $outFile);
print OUT "\#Node\tHighestCore\tStrengthThreshold\tCoreSize\n";

foreach my $node (sort {$core{$a} <=> $core{$b}} keys %core) {
	print OUT $node . "\t" . $core{$node} . "\t" . &digits($thresh{$node}, 4) . "\t" . $coreSize{$core{$node}} . "\n";
}
close OUT;



#------------------------------------------------------
sub acc {
# Arg: $numberToBeCutToAccuracy(ifNeeded)
# Ret: $accurateNumber
	sprintf('%.' . $accurateDigits . 'f', $_[0])
}

sub digits {
# Arg: ($numberToBeCutToNDigits, $N)
# Ret: $numberWithNDigits
	sprintf('%.' . $_[1] . 'f', $_[0])
}


sub fallback {
# Arg: (\%fallbackNetwork, \%fallbackStrength, \%network, \%nodeStrength).
	my $net = shift @_;
	my $str = shift @_;
	my $network = shift @_;
	my $strength = shift @_;

	%$net = ();
	%$str = ();

	foreach my $key1 (keys %$network) {
		foreach my $key2 (keys %{ $network->{$key1} }) {
			$net->{$key1}{$key2} = $network->{$key1}{$key2};
		}
	}

	%$str = %$strength;

	1;
}


sub loadNetwork {
# Arg: (\%network, \%nodeStrength, $fileName, $key1Col, $key2Col, $valCol).
	my $network = shift @_;
	my $strength = shift @_;
	my $file = shift @_;
	my $key1Col = shift @_;
	my $key2Col = shift @_;
	my $valCol = shift @_;

	open (FILE, '<' . $file) || die "Could not open $file ." . "\n";
	while (<FILE>) {
		next if ( m/^\s*\#/);
		next if ( m/^\s*$/ );
		s/^\s+//;
		chomp;

		my @split = split(/\t/,$_);

		my $val = abs($split[$valCol]);
		$val = 1 if $ifKCore eq "yes"; # strength = degree if $ifKCore eq "yes"

		# Edge A-B = B-A
		$network->{$split[$key1Col]}{$split[$key2Col]} = $val;
		$network->{$split[$key2Col]}{$split[$key1Col]} = $val;

		$strength->{$split[$key1Col]} = 0 if not exists $strength->{$split[$key1Col]};
		$strength->{$split[$key1Col]} += $val;
		$strength->{$split[$key2Col]} = 0 if not exists $strength->{$split[$key2Col]};
		$strength->{$split[$key2Col]} += $val;
	}
	close FILE;
	1;
}

sub min2 {
# Arg: ($val1, $val2)
# Ret: min(val1, val2)
	if (&acc($_[0]) <= &acc($_[1])) {
		return $_[0];
	} else {
		return $_[1];
	}
}

sub se {
	print STDERR $_[0] . "\n";
}

sub sortedLinkWeights {
# Arg: (\%@listWithSortedWeights, \%network).
	my $list = shift @_;
	my $network = shift @_;

	my %networkCheck;
	my $iterator = 0;
	foreach my $key1 (keys %$network) {
		foreach my $key2 (keys %{ $network->{$key1} } ) {
			unless ((exists $networkCheck{$key1}{$key2}) or (exists $networkCheck{$key2}{$key1})) {
				$list->[$iterator][0] = $key1;
				$list->[$iterator][1] = $key2;
				$list->[$iterator][2] = $network->{$key1}{$key2};
				$iterator++;
				$networkCheck{$key1}{$key2} = 1;
			}
		}
	}
	@$list = sort { $a->[2] <=> $b->[2] } @$list;
	1;
}

sub sortedListFromHash {
# Arg: (\%@listWithSortedValues, \%nodeStrength, $fileName, $key1Col, $key2Col, $valCol).
	my $list = shift @_;
	my $hashToBeSorted = shift @_;

	my $iterator = 0;
	foreach my $key (sort {$hashToBeSorted->{$a} <=> $hashToBeSorted->{$b}} keys %$hashToBeSorted) {
		$list->{$iterator} = $hashToBeSorted->{$key} if ref($list) eq "HASH";
		$list->[$iterator] = $hashToBeSorted->{$key} if ref($list) eq "ARRAY";
		$iterator++;
	}
	1;
}

sub updateHashes {
# Arg: (\%network, \%nodeStrength, \%core, \%thresh, \%checkNodes, $coreNumber, $node, $checkNode, $threshold)
	my $network = shift @_;
	my $strength = shift @_;
	my $core = shift @_;
	my $thresh = shift @_;
	my $checkNodes = shift @_;
	my $coreNumber = shift @_;
	my $node = shift @_;
	my $checkNode = shift @_;
	my $threshold = shift @_;

	$strength->{$node} -= $network->{$node}{$checkNode};
	# If a neighbor of $checkNode is 0, it can be removed right away as it doesn't have any links != 0 to other nodes
	if (&acc($strength->{$node}) == 0) {
		delete $strength->{$node};
		$core->{$node} = $coreNumber;
		$thresh->{$node} = $threshold;
	}
	# If a neighbor of $checkNode drops below $threshold, it must be removed, and its neighbors updated, in the next iteration
	elsif (&acc($strength->{$node}) <= &acc($threshold)) {
		$checkNodes->{$node} = 1;
	}
	# $checkNode is deleted along with its connections to $node, and is given a highest core number. 
	delete $network->{$node}{$checkNode};
	delete $network->{$checkNode}{$node};
	delete $strength->{$checkNode};
	delete $checkNodes->{$checkNode};
	$core->{$checkNode} = $coreNumber;
	$thresh->{$checkNode} = $threshold;
	1;
}


sub updateScorePlusHashes {
# Arg: (\%network, \%nodeStrength, $node1, $node2, $linkValue)
	my $network = shift @_;
	my $strength = shift @_;
	my $node1 = shift @_;
	my $node2 = shift @_;
	my $linkValue = shift @_;

	$strength->{$node1} -= $linkValue;
	$strength->{$node2} -= $linkValue;

	delete $network->{$node1}{$node2};
	delete $network->{$node2}{$node1};
	1;
}

