#!/usr/local/bin/perl -w

$lines_of_junk = 4;		# Number of lines of junk at top of .cm file
$nt = 21;			# Number of nucleotides in each molecule
$fields_of_junk = 4;		# Number of fields of junk at beginning of line
$angles = 11;			# Number of angles in each line

# List of numbers of files to be averaged. Files must be named "refine.#.cm".

#@files = (3, 5, 8, 11, 14, 17, 18, 20, 22, 23, 24, 27);
#@files = (3, 5, 11, 14, 17, 20, 22, 23, 24);
@files = (8, 18, 27);

$junk = undef;			# -w bait

# Read files into $angles

for ($file = 0; $file <= $#files; $file++) {
    open(CM, "refine." . sprintf("%2.2d", $files[$file]) . ".cm");
    for ($i = 1; $i <= $lines_of_junk; $i++) {
	$junk = <CM>;
    }
    for ($i = 1; $i <= $nt; $i++) {
	$angles[$file][$i] = [
	    (split(' ', <CM>))[
		$fields_of_junk - 1 .. $fields_of_junk + $angles - 1
	    ]
	];
	for ($angle = 0; $angle <= $angles - 1; $angle++) {
#	    $angles[$file][$i][$angle] +=
#		($angles[$file][$i][$angle] < 0 ? 360 : 0);
	}
    }
    close CM;
}

# Calculate and print averages and standard deviations

for ($i = 1; $i <= $nt; $i++) {
    for ($angle = 0; $angle <= $angles - 1; $angle++) {

	# Calculate average

	$sum = 0;
	for ($file = 0; $file <= $#files; $file++) {
	    $sum += $angles[$file][$i][$angle];
	}
	$average = $sum / ($#files + 1);

	# Calculate standard deviation

	$ssd = 0;
	for ($file = 0; $file <= $#files; $file++) {
	    $ssd += ($angles[$file][$i][$angle] - $average) ** 2;
	}
	$sd = sqrt($ssd/$#files);

	printf("%4.1d+%3.1d", $average, $sd);
#	printf("%4.1d+%3.1d", $average - ($average > 180 ? 360 : 0), $sd);
	print "\t" unless $angle == $angles;
    }
    print "\n";
}