File Coverage

blib/lib/HackaMol/Roles/ReadZmatRole.pm
Criterion Covered Total %
statement 18 93 19.3
branch 0 8 0.0
condition n/a
subroutine 6 8 75.0
pod 1 1 100.0
total 25 110 22.7


line stmt bran cond sub pod time code
1             $HackaMol::Roles::ReadZmatRole::VERSION = '0.053';
2             # ABSTRACT: Read files with molecular information
3             use Moose::Role;
4 11     11   5774 use HackaMol::PeriodicTable qw(%KNOWN_NAMES _trim);
  11         26  
  11         83  
5 11     11   48389 use Math::Vector::Real;
  11         29  
  11         1260  
6 11     11   70 use Carp;
  11         26  
  11         484  
7 11     11   68 use List::MoreUtils qw(singleton);
  11         39  
  11         610  
8 11     11   80  
  11         26  
  11         122  
9             with qw(
10             HackaMol::Roles::NERFRole
11             );
12              
13              
14             #xyz file and generate list of Atom object
15             my $self = shift;
16             my $fh = shift;
17 0     0 1   # my $file = shift;
18 0           # my $fh = FileHandle->new("<$file") or croak "unable to open $file";
19              
20             my @atoms;
21             my ( $n, $t ) = ( 0, 0 );
22 0          
23 0           my @zmat = <$fh>;
24             @zmat = _substitute_variables(@zmat);
25 0            
26 0           # we have 5 types of extensions
27             # A. SYM 0 x y z
28             # B. SYM
29             # C. SYM i R
30             # D. SYM i R j Ang
31             # E. SYM i R j Ang k Tors
32             # we need to filter the indices (can't lose the location)
33              
34             #type A
35             my @iA = grep { $zmat[$_] =~ m/^\s*\w+\s+0(\s+-*\d*\.*\d*){3}/ } 0 .. $#zmat;
36             my @inA = singleton( 0 .. $#zmat, @iA );
37 0            
  0            
38 0           #type B
39             my @iB = grep { $zmat[$_] =~ m/^\s*\w+\s*$/ } @inA;
40              
41 0           #type C
  0            
42             my @iC = grep { $zmat[$_] =~ m/^\s*\w+(\s+\d+\s+\d*\.*\d*)\s*$/ } @inA;
43              
44 0           #type D
  0            
45             my @iD = grep { $zmat[$_] =~ m/^\s*\w+(\s+\d+\s+\d*\.*\d*){2}\s*$/ } @inA;
46              
47 0           #type E
  0            
48             my @iE = grep {
49             $zmat[$_] =~ m/^\s*\w+(\s+\d+\s+\d*\.*\d*){2}\s+\d+\s+-*\d*\.*\d*\s*$/
50             } @inA;
51 0            
  0            
52             my $diff = @zmat - (@iA+@iB+@iC+@iD+@iE); #scalar context
53            
54 0           if ($diff){
55             print "Lines in Z-matrix: ", scalar (@zmat), " Number of lines to be processed: ", scalar (@zmat) - $diff, "\n";
56 0 0         print "Lines missed: ", $diff, "\n";
57 0           print "\n\nHere is your Z-matrix:\n";
58 0           print $_ foreach @zmat;
59 0           print "Indices of lines to be processed: ", join("\n", @iA, @iB, @iC, @iD, @iE);
60 0           croak "\nThere is something funky with your zmatrix";
61 0           }
62 0            
63             foreach my $ia (@iA) {
64             my ( $sym, $iat1, @xyz ) = split( ' ', $zmat[$ia] );
65 0           $atoms[$ia] = HackaMol::Atom->new(
66 0           name => $sym.$ia,
67 0           symbol => $sym,
68             coords => [ V(@xyz) ]
69             );
70             }
71            
72             foreach my $ib (@iB) {
73             my $sym = $zmat[$ib];
74 0           my $a = $self->init;
75 0           $sym =~ s/^\s+|\s+$//;
76 0           $atoms[$ib] = HackaMol::Atom->new(
77 0           name => $sym.$ib,
78 0           symbol => $sym,
79             coords => [$a]
80             );
81             }
82              
83             # print Dump 'B', \@atoms;
84              
85             foreach my $ic (@iC) {
86             my ( $sym, $iat1, $R ) = split( ' ', $zmat[$ic] );
87 0           my $a = $atoms[ $iat1 - 1 ]->xyz;
88 0           my $b = $self->extend_a( $a, $R );
89 0           $atoms[$ic] = HackaMol::Atom->new(
90 0           name => $sym.$ic,
91 0           symbol => $sym,
92             coords => [$b]
93             );
94             }
95              
96             # print Dump 'C', \@atoms;
97              
98             foreach my $id (@iD) {
99             my ( $sym, $iat1, $R, $iat2, $ang ) = split( ' ', $zmat[$id] );
100 0           my $a = $atoms[ $iat1 - 1 ]->xyz;
101 0           my $b = $atoms[ $iat2 - 1 ]->xyz;
102 0           my $c = $self->extend_ab( $b, $a, $R, $ang );
103 0           $atoms[$id] = HackaMol::Atom->new(
104 0           name => $sym.$id,
105 0           symbol => _trim($sym),
106             coords => [$c]
107             );
108             }
109              
110             # print Dump 'D', \@atoms;
111              
112             foreach my $ie (@iE) {
113             my ( $sym, $iat1, $R, $iat2, $ang, $iat3, $tor ) =
114 0           split( ' ', $zmat[$ie] );
115 0           my $a = $atoms[ $iat1 - 1 ]->xyz;
116             my $b = $atoms[ $iat2 - 1 ]->xyz;
117 0           my $c = $atoms[ $iat3 - 1 ]->xyz;
118 0           my $d = $self->extend_abc( $c, $b, $a, $R, $ang, $tor );
119 0           $atoms[$ie] = HackaMol::Atom->new(
120 0           name => $sym.$ie,
121 0           symbol => _trim($sym),
122             coords => [$d]
123             );
124             }
125             $atoms[$_]->iatom($_) foreach ( 0 .. $#atoms );
126             return (\@atoms);
127 0            
128 0           }
129              
130             my @Zmat = @_;
131              
132             chomp @Zmat;
133 0     0      
134             my %bin;
135 0           my %var = map {
136             my ($key,$val) = map{ s/^\s+|\s+$//; $_ } split(/\s*=\s*/,$_);
137 0           $bin{$key}++;
138             $key => $val,
139 0           } grep {/=/} @Zmat;
  0            
  0            
140 0            
141 0           # check for double entry of variables
142 0           my @too_many = grep {$bin{$_}>1} keys(%bin);
  0            
143             if (@too_many) {
144             carp "ReadZMatRole> you have more than one entry for these variables: ". join("\n", @too_many);
145 0           }
  0            
146 0 0          
147 0           @Zmat = grep {!/(^\#)|=|(^\s*$)/} @Zmat;
148              
149             foreach my $line (@Zmat){
150 0           my @vals = split (/ /, $line);
  0            
151             next unless @vals > 2;
152 0           $line = join(' ', $vals[0], map{ exists($var{$_}) ? $var{$_} : $_ } @vals[1 .. $#vals] );
153 0           }
154 0 0         return (@Zmat);
155 0 0         }
  0            
156              
157 0           no Moose::Role;
158              
159             1;
160 11     11   19763  
  11         44  
  11         70  
161              
162             =pod
163              
164             =head1 NAME
165              
166             HackaMol::Roles::ReadZmatRole - Read files with molecular information
167              
168             =head1 VERSION
169              
170             version 0.053
171              
172             =head1 SYNOPSIS
173              
174             my @atoms = HackaMol->new
175             ->read_zmat_atoms("some.zmat");
176              
177             =head1 DESCRIPTION
178              
179             The HackaMol::Roles::ReadZmatRole provides read_zmat_atoms for the flexible reading of Z-matrix files.
180             It supports inline cartesian coordinates and variables as in the following example:
181              
182             N 0 -12.781 3.620 15.274
183              
184             C 0 -11.976 4.652 15.944
185              
186             C 0 -12.722 6.019 15.985
187              
188             O 0 -13.133 6.378 14.897
189              
190             C 2 CBCA 3 CBCAC 4 CBCACO
191              
192             C 5 CBCA 2 CBCAC 3 CG1CBCAC
193              
194             C 5 CBCA 2 CBCAC 3 CG2CBCAC
195              
196             CBCA = 1.54
197              
198             CBCAC = 113.4
199              
200             CBCACO = 71.85
201              
202             CG1CBCAC = 54.
203              
204             CG2CBCAC = 180.
205              
206             =head1 METHODS
207              
208             =head2 read_zmat_atoms
209              
210             One argument: the filename
211             Returns a list of HackaMol::Atom objects.
212              
213             =head1 SEE ALSO
214              
215             =over 4
216              
217             =item *
218              
219             L<HackaMol>
220              
221             =item *
222              
223             L<HackaMol::Atom>
224              
225             =item *
226              
227             L<HackaMol::Roles::MolReadRole>
228              
229             =item *
230              
231             L<Protein Data Bank|http://pdb.org>
232              
233             =back
234              
235             =head1 CONSUMES
236              
237             =over 4
238              
239             =item * L<HackaMol::Roles::NERFRole>
240              
241             =back
242              
243             =head1 AUTHOR
244              
245             Demian Riccardi <demianriccardi@gmail.com>
246              
247             =head1 COPYRIGHT AND LICENSE
248              
249             This software is copyright (c) 2017 by Demian Riccardi.
250              
251             This is free software; you can redistribute it and/or modify it under
252             the same terms as the Perl 5 programming language system itself.
253              
254             =cut