Subversion Repositories wimsdev

Rev

Rev 14376 | Blame | Compare with Previous | Last modification | View Log | RSS feed

  1. #!/bin/sh
  2.  
  3. exec perl <<'EOF'
  4.  
  5. use strict "vars";
  6. use warnings;
  7.  
  8. push (@ARGV,split(' ', $ENV{'wims_exec_parm'})) if ($ENV{'wims_exec_parm'});
  9. my ($file)=shift @ARGV;
  10.  
  11. my @valence1=('Li','Na','K','F','Cl','Br','I','At');
  12. my @valence2=('Be','Mg','Ca','Sr','Ba','O','Se','Te','Mn','Hg','Zn','Cu','S');
  13. my @valence3=('B','Al','N','P','As','Fe');
  14. my @valence4=('C','Si','Ge','Sn','A','Q');
  15. my @other=('H');
  16. my ($nbatom, $nbbond);
  17. my $cnt=1; my $nH=0;
  18. my %hash = ();
  19. my $hash = \%hash;
  20. my %bond=();
  21. my $bond = \%bond;
  22. my %atom = ();
  23. my $atom = \%atom;
  24. my $cntmax=10000;
  25. my @a=();
  26. open(IN, "$file");
  27. while(<IN>) {
  28.   @a=split(" ", $_);
  29.   if($cnt==4){
  30.     $nbatom=$a[0]; $nbbond=$a[1];$cntmax=4+$nbatom;
  31.   };
  32.   if($cnt>4 && $cnt <= $cntmax) {
  33.     $hash->{$cnt-4}{name}=$a[3];
  34.     if(!$atom->{$a[3]}) { $atom->{$a[3]}=1} else { $atom->{$a[3]} ++}
  35.     if ( grep { $_ eq $a[3] } @valence4 ){ $hash->{$cnt-4}{hydrogen}=4}
  36.     elsif ( grep { $_ eq $a[3] } @valence1 ){ $hash->{$cnt-4}{hydrogen}=1}
  37.     elsif ( grep { $_ eq $a[3] } @valence2 ){ $hash->{$cnt-4}{hydrogen}=2}
  38.     elsif ( grep { $_ eq $a[3] } @valence3 ){ $hash->{$cnt-4}{hydrogen}=3}
  39.     elsif ( grep { $_ eq $a[3] } @other) {}
  40.     else {
  41.       #die( "$a[3]")
  42.     }
  43.   }
  44.   my $tt;
  45.   if ($cnt > $cntmax && $cnt <= $cntmax + $nbbond) {
  46.     if (!($hash->{$a[0]}{name} =~ /H\b/) && !($hash->{$a[1]}{name} =~ /H\b/)) {
  47.       $hash->{$a[0]}{'hydrogen'}-= $a[2];
  48.       $hash->{$a[1]}{'hydrogen'}-= $a[2];
  49.       $bond->{$cnt}=$a[0] . ',' . $a[1];
  50.       $tt=$cnt<= $cntmax ? $cnt: $cnt-$cntmax;
  51.       $bond->{$a[0]}=!$bond->{$a[0]} ? $tt: $bond->{$a[0]}.','. $tt;
  52.       $bond->{$a[1]}=!$bond->{$a[1]} ? $tt: $bond->{$a[1]}.','. $tt;
  53.     }
  54.     if (!($hash->{$a[0]}{name} =~ /H\b/) && ($hash->{$a[1]}{name} =~ /H\b/)) {
  55.       $hash->{$a[0]}{'hydrogenfichier'}++;
  56.     }
  57.     if (!($hash->{$a[1]}{name} =~ /H\b/) && ($hash->{$a[0]}{name} =~ /H\b/)) {
  58.       $hash->{$a[1]}{'hydrogenfichier'}++;
  59.     }
  60.   }
  61.   if ($_ =~s /M  CHG//) {
  62.     @a=split(" ", $_);
  63.     my $cc=1;
  64.     for my $ion (@a) {
  65.       if ($cc%2==0) {
  66.         $hash->{$ion}{'hydrogen'} += $a[$cc];
  67.       }
  68.       $cc++;
  69.     }
  70.   };
  71.   $cnt ++;
  72. }
  73. close IN;
  74. my $nbhydrogen='';
  75. for my $at (sort {$a<=>$b} keys %hash){
  76.   $hash{$at}{'hydrogen'}=0 if (!(defined $hash{$at}{'hydrogen'}) || $hash{$at}{'hydrogen'}<0);
  77.   $nH = $nH+ $hash{$at}{'hydrogen'};
  78.   $nbhydrogen .= "$at,$hash->{$at}{name},$hash->{$at}{'hydrogen'}\n";
  79. }
  80. my $nbhydrogenfichier='';
  81. for my $at (sort {$a<=>$b} keys %hash){
  82.   if(!($hash->{$at}{name} =~ /H\b/)){
  83.     $hash{$at}{'hydrogenfichier'}=0 if (!(defined $hash{$at}{'hydrogenfichier'})
  84.       || $hash{$at}{'hydrogenfichier'}<0);
  85.     $nbhydrogenfichier .= "$at,$hash->{$at}{name},$hash->{$at}{'hydrogenfichier'}\n";
  86.   }
  87. }
  88.  
  89. $atom{'H'}=$nH;
  90. my $formula=$atom{'C'} > 0 ? 'C' : '' ;
  91. $formula .= $atom{'C'} if $atom{'C'} > 1;
  92. $formula .= 'H' if $atom{'H'} > 0;
  93. $formula .= $atom{'H'} if $atom{'H'} > 1;
  94.  
  95. for my $at (sort keys %atom){
  96.   next if (!(defined ($atom{$at})) || $atom{$at}==0 );
  97.   next if ($at=~'C\b' || $at=~'H\b');
  98.   $formula .= $atom{$at}==1? $at : $at . $atom{$at};
  99. }
  100. my $bondprint="";
  101. for my $bd (sort {$a<=>$b} keys %bond){
  102.   next if (!(defined ($bond{$bd})));
  103.   if($bd <= $nbatom){
  104.     $bondprint .= "atom $bd,$bond{$bd}\n";
  105.   } else {
  106.     my $bd1=$bd-$nbatom-4; $bondprint .= "bond $bd1,$bond{$bd}\n";
  107.   }
  108. }
  109. print "$formula,[$nbhydrogen],[$nbhydrogenfichier],[$bondprint]";
  110.  
  111. EOF
  112.