Лаборатория изучения механизмов органических реакций

join_IRC

#!/usr/bin/perl -ws

#use Data::Dump 'pp';
#$Data::Dump::LINEWIDTH = 80;

our ($h,$help);

if ($h || $help) {
  (my $program = $0) =~ s/^.*[\/\\]//;
  print <<HELP;
Usage: $program IRC1.xyz IRC2.xyz ... | mldn
Зависимости: perl, скрипты cfxyz, xyz, reflect, renum, align
Стыкует IRC-траектории. Каждую следующую траекторию перенумеровывает по 
последней точке предыдущей и при необходимости обращает ее конфигурацию.
Нужно самому позаботиться, чтобы структуры на стыке траекторий совпадали.
-RMS=0.2  если при стыковке RMS будет больше, то die
-ref  структура или 1-я точка траектории, к которой подстраивать
      Если не задано, то подстраивается к первой траектории из аргументов 
HELP
  exit;
}
die "Usage: $0 [-RMS] [-ref=0.xyz] IRC1.xyz IRC2.xyz ... | mldn\n" unless @ARGV;

our ($RMS,$ref);

my $num = qr/-?\d+(?:\.\d+)?/;

$RMS ||= 0.2;

if ($ref) {
  my @out = `(xyz -n=1 $ref; xyz -n=1 $ARGV[0]) | cfxyz -all -vector -chiral`;
  my $vector1 = $out[0];
  chomp $vector1;
  die "Invalid vector between $ref and $ARGV[0]\n" if $vector1 !~ /^\d+(?:,\d+)+$/;
  my $rms1 = $1 if $out[1] =~ /RMS\s+($num)/;
  @out = `(xyz -n=1 $ref; xyz -n=1 $ARGV[0] | reflect) | cfxyz -all -vector -chiral`;
  my $vector2 = $out[0];
  chomp $vector2;
  die "Invalid vector between $ref and $ARGV[0]\n" if $vector2 !~ /^\d+(?:,\d+)+$/;
  my $rms2 = $1 if $out[1] =~ /RMS\s+($num)/;
  if ($rms1 <= $rms2) {
    die "RMS $rms1 between $ref and $ARGV[0]\n" if $rms1 > $RMS;
    `renum $vector1 $ARGV[0] > curr_IRC.xyz`;
  }
  else {
    die "RMS $rms2 beetwen $ref and $ARGV[0]\n" if $rms2 > $RMS;
    `reflect $ARGV[0] | renum $vector2 > curr_IRC.xyz`;
  }
}
else {
  if (@ARGV == 1) {
    `cat $ARGV[0]`;
    exit;
  }
  `cp $ARGV[0] curr_IRC.xyz`;
}

for (my $i=1; $i<@ARGV; $i++) {
  my @out = `(xyz -n=-1 curr_IRC.xyz; xyz -n=1 $ARGV[$i]) | cfxyz -all -vector -chiral`;
  my $vector1 = $out[0];
  chomp $vector1;
  die "Invalid vector between $ARGV[$i-1] and $ARGV[$i]\n" if $vector1 !~ /^\d+(?:,\d+)+$/;
  my $rms1 = $1 if $out[1] =~ /RMS\s+($num)/;
  @out = `(xyz -n=-1 curr_IRC.xyz; xyz -n=1 $ARGV[$i] | reflect) | cfxyz -all -vector -chiral`;
  my $vector2 = $out[0];
  chomp $vector2;
  die "Invalid vector between $ARGV[$i-1] and $ARGV[$i]\n" if $vector2 !~ /^\d+(?:,\d+)+$/;
  my $rms2 = $1 if $out[1] =~ /RMS\s+($num)/;
  if ($rms1 <= $rms2) {
    die "RMS $rms1 between $ARGV[$i-1] and $ARGV[$i]\n" if $rms1 > $RMS;
    `renum $vector1 $ARGV[$i] >> curr_IRC.xyz`;
  }
  else {
    die "RMS $rms2 beetwen $ARGV[$i-1] and $ARGV[$i]\n" if $rms2 > $RMS;
    `reflect $ARGV[$i] | renum $vector2 >> curr_IRC.xyz`;
  }
}

system "align curr_IRC.xyz";
unlink 'curr_IRC.xyz';