#!/usr/local/bin/perl
# 04.11.92
# Martin Schuetz, Institute of physical chemistry,
# University of Berne, Switzerland
# extracts geometry from G90 or CADPAC and produces SCHAKAL input
# file

die "You did not give G90 or CADPAC output file name as argument\n" if $#ARGV <
0;
die "You need only one argument: filename\n" if $#ARGV > 0;

$GCOUT = $ARGV[0];
$GCOUT =~ /([\w\/+#-]+)/;
$scal_fac = 1.0;
if ($#ARGV == 1) {
  $scal_fac = $ARGV[1];
  }
$OUT_root = $1;

@at_symbols = ('H ', 'He',
               'Li', 'Be', 'B ', 'C ', 'N ', 'O ', 'F ', 'Ne',
               'Na', 'Mg', 'Al', 'Si', 'P ', 'S ', 'Cl', 'Ar',
               'K ', 'Ca',
               'Ti', 'V ', 'Cr', 'Mn', 'Fe', 'Co', 'Ni', 'Cu', 'Zn',
                           'Ga', 'Ge', 'As', 'Se', 'Br', 'Kr', 'X ');

@at_weights = ( 1, 4,
                7, 9, 11, 12, 14, 16, 19, 20,
                23, 24, 27, 28, 31, 32, 35, 40 );

@ab_initio = ('unknown', 'Gaussian', 'CADPAC');

open(GCOUT,"<$GCOUT") || die "Could not open $GCOUT\n";

# check if Gaussian or CADPAC
$outf_type = do CHECK();
die "$GCOUT is neither Gaussian90 nor CADPAC outputfile" if !$outf_type;
printf "Hmm, this is a @ab_initio[$outf_type] output file\n";

# process output file...
do do_it_babe();
close(GCOUT);
printf "all done\n";

#
#
# -- subroutines
#
#

sub CHECK {
  local($chk) = 0;
  while (<GCOUT>) {
    if (/Gaussian/) { $chk = 1; }
    if (/CAMBRIDGE ANALYTIC DERIVATIVES PACKAGE/) { $chk = 2; }
    last if ($chk != 0);
    }
  $chk;
  }


sub do_it_babe {
  # go through output file, searching for "Nuclear coordinates"
  $out_no = 0;
  if ($outf_type == 1) { $search_exp = "Standard orientation:"; }
  if ($outf_type == 2) { $search_exp = "Nuclear coordinates"; }
  while (<GCOUT>) {
    if (/$search_exp/) {
      # skip next 4 lines
      for ($i = 0; $i < 4; $i++ ) {
        $line = <GCOUT>;
      }
      $out_no++;
      # open output file
      $temp = sprintf("%d",$out_no);
      $OUTF = $OUT_root . '.geo' . $temp;
      open(OUTF,">$OUTF") || die "Could not open $OUTF\n";
      printf(OUTF "CELL\n");
      $n_at = 0;
      if ($outf_type == 1) { $serch_exp = "-------------"; }
      if ($outf_type == 2) { $serch_exp = "BOND LENGTHS IN BOHR"; }          
      while ($line = <GCOUT>) {
        last if ($line =~ /$serch_exp/);  # when last line
        $n_at++;
        $line =~ /\S+\s+(\S+)\s+(\S+)\s+(\S+)\s+(\S+)/;
        $x_coor = $2;
        $y_coor = $3;
        $z_coor = $4;
        if ($outf_type == 1) {
          if($1 < $#at_symbols) {
            $at_symb = $at_symbols[$1-1];
            }
          else {
            $at_symb = $at_symbols[$#at_symbols];
            }
          }
        if ($outf_type == 2) { 
          $at_symb = $1;
          $x_coor *= 0.529177249;
          $y_coor *= 0.529177249;
          $z_coor *= 0.529177249;
          }
        printf(OUTF "ATOM  %5s  %9.4f %9.4f %9.4f\n", $at_symb, 
               $x_coor, $y_coor, $z_coor);
        }
      close(OUTF);
      }
    }
  }