+ if (($1 < 1) || ($1 > $data->{ndihedraltypes}));
+ ++$i;
+ $data->{middlebondtorsioncoeff}[$j] = $2;
+ next;
+ }
+
+ die "Too many entries in MiddleBondTorsion Coeffs section"
+ if ($i > $data->{ndihedraltypes});
+ die "Too few entries in MiddleBondTorsion Coeffs section"
+ if ($i < $data->{ndihedraltypes});
+ die "Multiple dihedral coefficient assignments to the same dihedral type"
+ if (scalar @{$data->{middlebondtorsioncoeff}} != $data->{ndihedraltypes});
+
+ last;
+ }
+ } elsif ($1 eq "BondBond13 Coeffs") {
+ $data->{bondbond13coeff} = [];
+ $i = 0;
+ section_check($fh,$data,"BondBond13 Coeffs");
+
+ while (get_next($fh)) {
+
+ if (/^\s*([0-9]+)\s+([-+.eE0-9 ]+)\s*$/) {
+ $j = $1 - 1;
+ die "BondBond13 type $1 is out of range"
+ if (($1 < 1) || ($1 > $data->{ndihedraltypes}));
+ ++$i;
+ $data->{bondbond13coeff}[$j] = $2;
+ next;
+ }
+
+ die "Too many entries in BondBond13 Coeffs section"
+ if ($i > $data->{ndihedraltypes});
+ die "Too few entries in BondBond13 Coeffs section"
+ if ($i < $data->{ndihedraltypes});
+ die "Multiple dihedral coefficient assignments to the same dihedral type"
+ if (scalar @{$data->{bondbond13coeff}} != $data->{ndihedraltypes});
+
+ last;
+ }
+ } elsif ($1 eq "AngleTorsion Coeffs") {
+ $data->{angletorsioncoeff} = [];
+ $i = 0;
+ section_check($fh,$data,"AngleTorsion Coeffs");
+
+ while (get_next($fh)) {
+
+ if (/^\s*([0-9]+)\s+([-+.eE0-9 ]+)\s*$/) {
+ $j = $1 - 1;
+ die "AngleTorsion type $1 is out of range"
+ if (($1 < 1) || ($1 > $data->{ndihedraltypes}));
+ ++$i;
+ $data->{angletorsioncoeff}[$j] = $2;
+ next;
+ }
+
+ die "Too many entries in AngleTorsion Coeffs section"
+ if ($i > $data->{ndihedraltypes});
+ die "Too few entries in AngleTorsion Coeffs section"
+ if ($i < $data->{ndihedraltypes});
+ die "Multiple dihedral coefficient assignments to the same dihedral type"
+ if (scalar @{$data->{angletorsioncoeff}} != $data->{ndihedraltypes});
+
+ last;
+ }
+ } elsif ($1 eq "Improper Coeffs") {
+ $data->{impropercoeff} = [];
+ $i = 0;
+ section_check($fh,$data,"Improper Coeffs");
+
+ while (get_next($fh)) {
+
+ if (/^\s*([0-9]+)\s+([-+.eE0-9 ]+)\s*$/) {
+ $j = $1 - 1;
+ die "Improper type $1 is out of range"
+ if (($1 < 1) || ($1 > $data->{nimpropertypes}));
+ ++$i;
+ $data->{impropercoeff}[$j] = $2;
+ next;
+ }
+
+ die "Too many entries in Improper Coeffs section"
+ if ($i > $data->{nimpropertypes});
+ die "Too few entries in Improper Coeffs section"
+ if ($i < $data->{nimpropertypes});
+ die "Multiple improper coefficient assignments to the same improper type"
+ if (scalar @{$data->{impropercoeff}} != $data->{nimpropertypes});
+
+ last;
+ }
+ } elsif ($1 eq "AngleAngle Coeffs") {
+ $data->{angleanglecoeff} = [];
+ $i = 0;
+ section_check($fh,$data,"AngleAngle Coeffs");
+
+ while (get_next($fh)) {
+
+ if (/^\s*([0-9]+)\s+([-+.eE0-9 ]+)\s*$/) {
+ $j = $1 - 1;
+ die "AngleAngle type $1 is out of range"
+ if (($1 < 1) || ($1 > $data->{nimpropertypes}));
+ ++$i;
+ $data->{angleanglecoeff}[$j] = $2;
+ next;
+ }
+
+ die "Too many entries in AngleAngle Coeffs section"
+ if ($i > $data->{nimpropertypes});
+ die "Too few entries in AngleAngle Coeffs section"
+ if ($i < $data->{nimpropertypes});
+ die "Multiple angleangle coefficient assignments to the same angle type"
+ if (scalar @{$data->{angleanglecoeff}} != $data->{nimpropertypes});
+
+ last;
+ }
+ } elsif ($1 eq "Atoms") {
+ $data->{tag} = [];
+ $data->{type} = [];
+ $data->{molid} = [];
+ $data->{charge} = [];
+ $data->{posx} = [];
+ $data->{posy} = [];
+ $data->{posz} = [];
+ $data->{imgx} = [];
+ $data->{imgy} = [];
+ $data->{imgz} = [];
+ $i = 0;
+ section_check($fh,$data,"Atoms");
+
+ while (get_next($fh)) {
+ if (/^\s*([0-9]+)\s+([0-9]+)\s+([0-9]+)\s+([-+.eE0-9]+)\s+([-+.eE0-9]+)\s+([-+.eE0-9]+)\s+([-+.eE0-9]+)(|\s+(-?[0-9]+)\s+(-?[0-9]+)\s+(-?[0-9]+))\s*$/) {
+
+ $k = $1 - 1;
+ die "Atom id $1 is out of range"
+ if (($1 < 1) || ($1 > $data->{natoms}));
+
+ $j = $3 - 1;
+ die "Atom type $2 is out of range"
+ if (($3 < 1) || ($3 > $data->{natomtypes}));
+
+ ++$i;
+ $data->{tag}[$k] = $1;
+ $data->{molid}[$k] = $2;
+ $data->{type}[$k] = $3;
+ $data->{charge}[$k] = $4;
+ $data->{posx}[$k] = $5;
+ $data->{posy}[$k] = $6;
+ $data->{posz}[$k] = $7;
+ $data->{imgx}[$k] = 0;
+ $data->{imgy}[$k] = 0;
+ $data->{imgz}[$k] = 0;
+ if (! $8 eq "") {
+ $data->{imgx}[$k] = $9;
+ $data->{imgy}[$k] = $10;
+ $data->{imgz}[$k] = $11;
+ }
+ next;
+# } else {
+# print "Atoms: $_\n";
+ }
+
+ die "Too many entries in Atoms section: $i vs. $data->{natoms}"
+ if ($i > $data->{natoms});
+ die "Too few entries in Atoms section: $i vs. $data->{natoms}"
+ if ($i < $data->{natoms});
+ die "Multiple atoms assigned to the same atom ID"
+ if (scalar @{$data->{tag}} != $data->{natoms});
+
+ last;
+ }
+ } elsif ($1 eq "Velocities") {
+ $data->{velx} = [];
+ $data->{vely} = [];
+ $data->{velz} = [];
+ $i = 0;
+ section_check($fh,$data,"Velocities");
+
+ while (get_next($fh)) {
+ if (/^\s*([0-9]+)\s+([-+.eE0-9]+)\s+([-+.eE0-9]+)\s+([-+.eE0-9]+)\s*$/) {
+
+ $k = $1 - 1;
+ die "Atom id $1 is out of range"
+ if (($1 < 1) || ($1 > $data->{natoms}));
+
+ ++$i;
+ $data->{velx}[$k] = $2;
+ $data->{vely}[$k] = $3;
+ $data->{velz}[$k] = $4;
+ next;
+ }
+
+ die "Too many entries in Velocities section"
+ if ($i > $data->{natoms});
+ die "Too few entries in Velocities section"
+ if ($i < $data->{natoms});
+ die "Multiple velocities assigned to the same atom ID"
+ if (scalar @{$data->{velx}} != $data->{natoms});
+
+ last;
+ }
+ } elsif ($1 eq "Bonds") {
+ $data->{bondt} = [];
+ $data->{bond1} = [];
+ $data->{bond2} = [];
+ $i = 0;
+ section_check($fh,$data,"Bonds");
+
+ while (get_next($fh)) {
+ if (/^\s*([0-9]+)\s+(-?[0-9]+)\s+([0-9]+)\s+([0-9]+)\s*$/) {
+
+ $k = $1 - 1;
+ die "Bond id $1 is out of range"
+ if (($1 < 1) || ($1 > $data->{nbonds}));
+
+ die "Bond type $2 is out of range"
+ if (($2 == 0) || ($2 > $data->{nbondtypes}));
+
+ die "Bond atom 1 ID $3 is out of range"
+ if (($3 < 1) || ($3 > $data->{natoms}));
+ die "Bond atom 2 ID $4 is out of range"
+ if (($4 < 1) || ($4 > $data->{natoms}));
+
+ ++$i;
+ $data->{bondt}[$k] = $2;
+ $data->{bond1}[$k] = $3;
+ $data->{bond2}[$k] = $4;
+ next;
+ }
+
+ die "Too many entries in Bonds section"
+ if ($i > $data->{nbonds});
+ die "Too few entries in Bonds section"
+ if ($i < $data->{nbonds});
+ die "Multiple bonds assigned to the same bond ID"
+ if (scalar @{$data->{bondt}} != $data->{nbonds});
+
+ last;
+ }
+ } elsif ($1 eq "Angles") {
+ $data->{anglet} = [];
+ $data->{angle1} = [];
+ $data->{angle2} = [];
+ $data->{angle3} = [];
+ $i = 0;
+ section_check($fh,$data,"Angles");
+
+ while (get_next($fh)) {
+ if (/^\s*([0-9]+)\s+(-?[0-9]+)\s+([0-9]+)\s+([0-9]+)\s+([0-9]+)\s*$/) {
+
+ $k = $1 - 1;
+ die "Angle id $1 is out of range"
+ if (($1 < 1) || ($1 > $data->{nangles}));
+
+ die "Angle type $2 is out of range"
+ if (($2 == 0) || ($2 > $data->{nangletypes}));
+
+ die "Angle atom 1 ID $3 is out of range"
+ if (($3 < 1) || ($3 > $data->{natoms}));
+ die "Angle atom 2 ID $4 is out of range"
+ if (($4 < 1) || ($4 > $data->{natoms}));
+ die "Angle atom 3 ID $5 is out of range"
+ if (($5 < 1) || ($5 > $data->{natoms}));
+
+ ++$i;
+ $data->{anglet}[$k] = $2;
+ $data->{angle1}[$k] = $3;
+ $data->{angle2}[$k] = $4;
+ $data->{angle3}[$k] = $5;
+ next;
+ }
+
+ die "Too many entries in Angles section"
+ if ($i > $data->{nangles});
+ die "Too few entries in Angles section"
+ if ($i < $data->{nangles});
+ die "Multiple angles assigned to the same angle ID"
+ if (scalar @{$data->{anglet}} != $data->{nangles});
+
+ last;
+ }
+ } elsif ($1 eq "Dihedrals") {
+ $data->{dihedralt} = [];
+ $data->{dihedral1} = [];
+ $data->{dihedral2} = [];
+ $data->{dihedral3} = [];
+ $data->{dihedral4} = [];
+ $i = 0;
+ section_check($fh,$data,"Dihedrals");
+
+ while (get_next($fh)) {
+ if (/^\s*([0-9]+)\s+(-?[0-9]+)\s+([0-9]+)\s+([0-9]+)\s+([0-9]+)\s+([0-9]+)\s*$/) {
+
+ $k = $1 - 1;
+ die "Dihedral id $1 is out of range"
+ if (($1 < 1) || ($1 > $data->{ndihedrals}));
+
+ die "Dihedral type $2 is out of range"
+ if (($2 == 0) || ($2 > $data->{ndihedraltypes}));
+
+ die "Dihedral atom 1 ID $3 is out of range"
+ if (($3 < 1) || ($3 > $data->{natoms}));
+ die "Dihedral atom 2 ID $4 is out of range"
+ if (($4 < 1) || ($4 > $data->{natoms}));
+ die "Dihedral atom 3 ID $5 is out of range"
+ if (($5 < 1) || ($5 > $data->{natoms}));
+ die "Dihedral atom 4 ID $6 is out of range"
+ if (($6 < 1) || ($6 > $data->{natoms}));
+
+ ++$i;
+ $data->{dihedralt}[$k] = $2;
+ $data->{dihedral1}[$k] = $3;
+ $data->{dihedral2}[$k] = $4;
+ $data->{dihedral3}[$k] = $5;
+ $data->{dihedral4}[$k] = $6;
+ next;
+ }
+
+ die "Too many entries in Dihedrals section"
+ if ($i > $data->{ndihedrals});
+ die "Too few entries in Dihedrals section"
+ if ($i < $data->{ndihedrals});
+ die "Multiple dihedrals assigned to the same dihedral ID"
+ if (scalar @{$data->{dihedralt}} != $data->{ndihedrals});
+
+ last;
+ }
+ } elsif ($1 eq "Impropers") {
+ $data->{impropert} = [];
+ $data->{improper1} = [];
+ $data->{improper2} = [];
+ $data->{improper3} = [];
+ $data->{improper4} = [];
+ $i = 0;
+ section_check($fh,$data,"Impropers");
+
+ while (get_next($fh)) {
+ if (/^\s*([0-9]+)\s+(-?[0-9]+)\s+([0-9]+)\s+([0-9]+)\s+([0-9]+)\s+([0-9]+)\s*$/) {
+
+ $k = $1 - 1;
+ die "Improper id $1 is out of range"
+ if (($1 < 1) || ($1 > $data->{nimpropers}));
+
+ die "Improper type $2 is out of range"
+ if (($2 == 0) || ($2 > $data->{nimpropertypes}));
+
+ die "Improper atom 1 ID $3 is out of range"
+ if (($3 < 1) || ($3 > $data->{natoms}));
+ die "Improper atom 2 ID $4 is out of range"
+ if (($4 < 1) || ($4 > $data->{natoms}));
+ die "Improper atom 3 ID $5 is out of range"
+ if (($5 < 1) || ($5 > $data->{natoms}));
+ die "Improper atom 4 ID $6 is out of range"
+ if (($6 < 1) || ($6 > $data->{natoms}));
+
+ ++$i;
+ $data->{impropert}[$k] = $2;
+ $data->{improper1}[$k] = $3;
+ $data->{improper2}[$k] = $4;
+ $data->{improper3}[$k] = $5;
+ $data->{improper4}[$k] = $6;
+ next;
+ }
+
+ die "Too many entries in Impropers section"
+ if ($i > $data->{nimpropers});
+ die "Too few entries in Impropers section"
+ if ($i < $data->{nimpropers});
+ die "Multiple impropers assigned to the same improper ID"
+ if (scalar @{$data->{impropert}} != $data->{nimpropers});
+
+ last;
+ }
+ } else {
+ die "Bad data: $_";
+ }
+
+ last unless ($_);
+ } else {
+ die "Bad data: $_";
+ }
+ }
+}
+
+sub floatdiff {
+ my ($n1,$n2,$rel) = @_;
+
+ my $diff = abs($n1-$n2);
+ my $avg = (abs($n1)+abs($n2))*0.5;
+ return 0 if ($avg == 0.0);
+ if ($rel) {
+# print "relative difference: ",$diff/$avg," vs. $small\n";
+ return 0 if ($diff/$avg < $small);
+ } else {
+# print "absolute difference: ",$diff," vs. $small\n";