#!/usr/bin/perl

# Extract all the energies from a minimization run, all to separate files.
# Allows for an optional prefix, otherwise they go to summary.VARNAME.
# Usage: process-minimization-mdout.pl [mdout file] [output prefix]

# List of things we want extracted and what storage variable we'll
#  use while iterating.
%variable_mapping = ('ENERGY' => 'energy',
		     'RMS' => 'rms',
		     'BOND' => 'bond',
		     'ANGLE' => 'angle',
		     'DIHED' => 'dihedral',
		     'VDWAALS' => 'vdwaals',
		     'EEL' => 'eel',
		     'HBOND' => 'hbond',
		     'VDW14' => 'vdw14',
		     'EEL14', => 'eel14',
		     'RESTRAINT' => 'restraint');

$ARGC = @ARGV;

unless ($ARGC >= 1) {die("Usage: $0 <mdout file> [output prefix]\n");}

# Try to set output prefix.
if ($ARGV[1]) {
  $prefix = $ARGV[1];
} else {
  $prefix = $ARGV[0];
}

# Attempt input.
open (INPUT, "<$ARGV[0]") || die("Can't open input $ARGV[0]: $!\n");

# Iterate over input until a meaningful line is reached.
# As we find each energy summary, cache it in an appropriate hash.
while (<INPUT>) {

  # Detect the start of a block.
  if (/NSTEP.*ENERGY.*RMS.*GMAX/) {
    # Get next line.
    $_ = <INPUT>;
    
    # Process run info.
    ($step, $energy, $rms) = /\s+(\d+)\s+(.\d+\.\d+E[-+]\d+)\s+(.\d+\.\d+E[-+]\d+)\s+/;

    # Eat two lines and start block-processing.
    $_ = <INPUT>; $_ = <INPUT>;

    if (/BOND.*ANGLE.*DIHED/) {
      ($bond, $angle, $dihedral) =
	/BOND.*=\s*(.*\d*\.\d*).*ANGLE.*=\s*(.*\d*\.\d*).*DIHED.*=\s*(.*\d*\.\d*)/;
      $_ = <INPUT>;
    }

    if (/VDWAALS.*EEL.*HBOND/) {
      ($vdwaals,$eel,$hbond) =
	/VDWAALS.*=\s*(.*\d*\.\d*).*EEL.*=\s*(.*\d*\.\d*).*HBOND.*=\s*(.*\d*\.\d*)/;

      $_ = <INPUT>;
    }

    if (/1-4 VDW.*1-4 EEL.*RESTRAINT/) {
      ($vdw14,$eel14,$restraint) =
	/1-4 VDW.*=\s*(.*\d*\.\d*).*1-4 EEL.*=\s*(.*\d*\.\d*).*RESTRAINT.*=\s*(.*\d*\.\d*)/;

      $_ = <INPUT>;
    }

    # Now that block processing is finished, update all arrays.
    # Note use of clever referencing trick.
    foreach $arrayname (keys(%variable_mapping)) {
      $$arrayname{$step} = eval "\$$variable_mapping{$arrayname}";
      
    }
  }

  # Continue to process blocks.
}


# Now that one has the energy fully read, spew it.
foreach $arrayname (keys(%variable_mapping)) {

  open (OUT, ">$prefix.$arrayname"); # Fail silently.

  %outarray = eval "\%$arrayname";

  foreach $time (sort {$a <=> $b} keys(%outarray)) {
    print OUT "$time $outarray{$time}\n";
  }

  close(OUT);
}

exit;