From 6e719f2d94699d8c6d1533766db18f95d71496e5 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Mon, 3 Oct 2016 07:07:28 -0400 Subject: [PATCH] remove trailing whitespace --- tools/ch2lmp/README | 30 ++--- tools/ch2lmp/charmm2lammps.pl | 230 +++++++++++++++++----------------- tools/ch2lmp/lammps2pdb.pl | 194 ++++++++++++++-------------- 3 files changed, 227 insertions(+), 227 deletions(-) diff --git a/tools/ch2lmp/README b/tools/ch2lmp/README index fa031bd757..25c0a7a9a0 100644 --- a/tools/ch2lmp/README +++ b/tools/ch2lmp/README @@ -41,9 +41,9 @@ should be -option=value (e.g. -border=5). -border add border to all sides of simulation box [default: 0 A] -ax rotation around x-axis -ay rotation around y-axis --az rotation around z-axis +-az rotation around z-axis -cd correction for dihedral for carbohydrate systems --cmap add CMAP section to data file and fix cmap command lines in +-cmap add CMAP section to data file and fix cmap command lines in input script" (NOTE: requires use of *.pdb file) In the "example" folder, you will find example files that were created @@ -69,26 +69,26 @@ file has to the corresponding names in the charmm FF files. You'll need to add a "pdbalias residue x xnew" line for each change that needs to be made. The *.pgn should contain something like this: -package require psfgen +package require psfgen topology top_all27_na.rtf -pdbalias residue A ADE -pdbalias residue T THY -pdbalias residue G GUA -pdbalias residue C CYT +pdbalias residue A ADE +pdbalias residue T THY +pdbalias residue G GUA +pdbalias residue C CYT . . . -segment A {pdb 1ac7_pared.pdb} -coordpdb 1ac7_pared.pdb A -guesscoord -writepdb 1ac7.pdb -writepsf charmm 1ac7.psf -exit +segment A {pdb 1ac7_pared.pdb} +coordpdb 1ac7_pared.pdb A +guesscoord +writepdb 1ac7.pdb +writepsf charmm 1ac7.psf +exit 5) Type "vmd -e 1ac7.pgn" to build the 1ac7.psf file, and the new 1ac7.pdb file. -6) Run charmm2lammps.pl by typing: +6) Run charmm2lammps.pl by typing: "perl charmm2lammps.pl all27_na 1ac7 -charmm -border=1 -pdb_ctrl -water -ions" 7) Run lammps by typing: "lmp < 1ac7.in" @@ -105,7 +105,7 @@ molecule. The -pdb_ctrl option produces the 1ac7_ctrl.pdb file that can be visualized in a standard visualization package such as VMD. The -charmm option put comments into the LAMMPS data file (everything after the # sign is a comment) for user convenience in tracking atom -types etc. according to CHARMM nomenclature. +types etc. according to CHARMM nomenclature. The example molecule provided above (i.e., 1ac7) is a DNA fragment. If instead, a peptide longer than 2 amino acid residues or a protein diff --git a/tools/ch2lmp/charmm2lammps.pl b/tools/ch2lmp/charmm2lammps.pl index d02b68eae7..dbc6ba003c 100755 --- a/tools/ch2lmp/charmm2lammps.pl +++ b/tools/ch2lmp/charmm2lammps.pl @@ -10,7 +10,7 @@ # 20050212 Needed (in the same directory): # - $project.crd ; Assumed to be correct and running # - $project.psf ; CHARMM configs -# - top_$forcefield.rtf ; +# - top_$forcefield.rtf ; # - par_$forcefield.prm ; # Ouput: # - $project.data ; LAMMPS data file @@ -41,10 +41,10 @@ # # General Many thanks to Paul S. Crozier for checking script validity # against his projects. -# Also thanks to Xiaohu Hu (hux2@ornl.gov) and Robert A. Latour -# (latourr@clemson.edu), David Hyde-Volpe, and Tigran Abramyan, -# Clemson University and Chris Lorenz (chris.lorenz@kcl.ac.uk), -# King's College London for their efforts to add CMAP sections, +# Also thanks to Xiaohu Hu (hux2@ornl.gov) and Robert A. Latour +# (latourr@clemson.edu), David Hyde-Volpe, and Tigran Abramyan, +# Clemson University and Chris Lorenz (chris.lorenz@kcl.ac.uk), +# King's College London for their efforts to add CMAP sections, # which is implemented using the option flag "-cmap". # Initialization @@ -52,7 +52,7 @@ sub Test { my $name = shift(@_); - + printf("Error: file %s not found\n", $name) if (!scalar(stat($name))); return !scalar(stat($name)); } @@ -63,7 +63,7 @@ my $k = 0; my @dir = ("x", "y", "z"); my @options = ("-help", "-nohints", "-water", "-ions", "-center", - "-quiet", "-pdb_ctrl", "-l", "-lx", "-ly", "-lz", + "-quiet", "-pdb_ctrl", "-l", "-lx", "-ly", "-lz", "-border", "-ax", "-ay", "-az", "-cmap"); my @remarks = ("display this message", "do not print type and style hints in data file", @@ -83,7 +83,7 @@ "generate a CMAP section in data file" ); my $notes; - + $program = "charmm2lammps"; $version = "1.9.0"; $year = "2016"; @@ -99,7 +99,7 @@ $L = (0, 0, 0); $cmap = 0; @R = M_Unit(); - + $notes = " * The average of extremes is used as the origin\n"; $notes .= " * Residues are numbered sequentially\n"; $notes .= " * Water is added on an FCC lattice: allow 5 ps for"; @@ -119,7 +119,7 @@ $notes .= " - project.data LAMMPS data file\n"; $notes .= " - project.in suggested LAMMPS input script\n"; $notes .= " - project_ctrl.pdb control file when requested\n"; - + foreach (@ARGV) { if (substr($_, 0, 1) eq "-") @@ -200,7 +200,7 @@ return "{".$v[0].", ".$v[1].", ".$v[2]."}"; } - + sub V_Add { my @v1 = splice(@_, 0, 3); @@ -208,8 +208,8 @@ return ($v1[0]+$v2[0], $v1[1]+$v2[1], $v1[2]+$v2[2]); } - - + + sub V_Subtr { my @v1 = splice(@_, 0, 3); @@ -217,8 +217,8 @@ return ($v1[0]-$v2[0], $v1[1]-$v2[1], $v1[2]-$v2[2]); } - - + + sub V_Dot { my @v1 = splice(@_, 0, 3); @@ -226,8 +226,8 @@ return $v1[0]*$v2[0]+$v1[1]*$v2[1]+$v1[2]*$v2[2]; } - - + + sub V_Mult { my @v = splice(@_, 0, 3); @@ -236,12 +236,12 @@ return ($f*$v[0], $f*$v[1], $f*$v[2]); } - + sub M_String { my $string; - - for (my $i=0; $i<3; ++$i) + + for (my $i=0; $i<3; ++$i) { $string .= ", " if ($i); $string .= V_String(splice(@_, 0, 3)); @@ -258,7 +258,7 @@ @_[2], @_[5], @_[8]); } - + sub M_Dot { my @v11 = splice(@_, 0, 3); @@ -268,7 +268,7 @@ my @v21 = splice(@m, 0, 3); my @v22 = splice(@m, 0, 3); my @v23 = splice(@m, 0, 3); - + return ( V_Dot(@v11, @v21), V_Dot(@v11, @v22), V_Dot(@v11, @v23), V_Dot(@v12, @v21), V_Dot(@v12, @v22), V_Dot(@v12, @v23), @@ -279,7 +279,7 @@ sub M_Unit { return (1,0,0, 0,1,0, 0,0,1); } sub PI { return 4*atan2(1,1); } - + sub M_Rotate { # vmd convention my $n = shift(@_); @@ -294,7 +294,7 @@ return ($cos,-$sin,0, $sin,$cos,0, 0,0,1) if ($n==2); # around z-axis return M_Unit(); } - + sub MV_Dot { @@ -305,10 +305,10 @@ return (V_Dot(@v11, @v2), V_Dot(@v12, @v2), V_Dot(@v13, @v2)); } - - + + # CHARMM input - + sub PSFConnectivity { my $n = PSFGoto(bonds); @@ -353,8 +353,8 @@ $dihedral_flag = 1; return $ndihedral; } - - + + sub CreatePSFIndex # make an index of id { # locations my @psf_ids = ("!NATOM","!NBOND:","!NTHETA:","!NPHI:","!NIMPHI:"); @@ -374,13 +374,13 @@ } } - + sub PSFGoto # goto $ident in { CreatePSFIndex() if (!scalar(%PSFIndex)); my $id = shift(@_); my @n = split(" ", $PSFIndex{$id}); - + @PSFBuffer = (); # return PSFDihedrals() if ($id eq "dihedrals"); if (!scalar(@n)) @@ -415,10 +415,10 @@ { my $items = shift(@_); my $n = $items; - + if ($psf_ncols>7) { printf(PSF_CTRL "\n"); $psf_ncols = 0; } - foreach(@_) - { + foreach(@_) + { printf(PSF_CTRL " %7d", $_); ++$psf_ncols; if ((!--$n) && ($psf_ncols>7)) @@ -434,7 +434,7 @@ sub CRDGoto { my $n; - + return if (shift(@_) ne "atoms"); open(CRD, "<".($pdb ? $Pdb : $Crd)) if (fileno(CRD) eq ""); seek(CRD, 0, SEEK_SET); @@ -451,20 +451,20 @@ my @data = (); my $c = 0; my $line; - + while (substr($line = , 0, 4) ne "ATOM") {}; chop($line); foreach (@n) { push(@data, substr($line, ($c += $_)-$_, $_)); } return @data[1, 8, 5, 3, 11, 12, 13, 17, 8, 15]; } - + sub Delete { my $item = shift(@_); my $k = 0; my @list; - + foreach (@_) { my @tmp = split(" "); @@ -483,9 +483,9 @@ my $flag = $list[0] gt $list[-1]; my $j = $n; my $tmp; - + return "" if (scalar(@list)<$n); - $flag = $list[1] gt $list[-2] + $flag = $list[1] gt $list[-2] if ((scalar(@list)>3)&&($list[0] eq $list[-1])); for (my $i=0; $i<$n; ++$i) { @@ -506,7 +506,7 @@ for (my $i=0; $i<$n; ++$i) { my @tmp = split(" ", ); - $tmp[5] = $symbols{$tmp[5]} + $tmp[5] = $symbols{$tmp[5]} if ((substr($tmp[5],0,1) lt '0')||(substr($tmp[5],0,1) gt '9')); push(@atom_types, $tmp[5]); ++$list{$tmp[5]}; @@ -524,7 +524,7 @@ return sort({$a<=>$b} keys(%list)); } - + sub Markers { my %markers = ( @@ -543,14 +543,14 @@ { my @cols = @_; my $f = (scalar(@cols)>3)&&(substr($cols[3],0,1) ne "!"); - my @tmp = (-$cols[1], $cols[2], + my @tmp = (-$cols[1], $cols[2], $f ? -$cols[4]:-$cols[1], $f ? $cols[5]:$cols[2]); $tmp[1] *= 2.0**(5/6); # adjust sigma $tmp[3] *= 2.0**(5/6); # adjust sigma 1-4 return join(" ", @tmp); } - + sub AtomParameters # non-bonded parameters { my @types; @@ -574,7 +574,7 @@ if ($markers{$cols[0]} ne "") { $read = ($markers{$cols[0]} eq "0") ? 1 : 0; } } - $list[$types{HT}] = NonBond(0, -0.046, 0.2245) + $list[$types{HT}] = NonBond(0, -0.046, 0.2245) if ($water_dens&&($list[$types{HT}] eq "")); $list[$types{OT}] = NonBond(0, -0.152100, 1.768200) if ($water_dens&&($list[$types{OT}] eq "")); @@ -593,7 +593,7 @@ my $id = (bonds, angles, dihedrals, impropers)[$mode]; my $n = PSFGoto($id); my %list; - + for (my $i=0; $i<$n; ++$i) { my @tmp = (); @@ -624,7 +624,7 @@ { # my $mode = shift(@_); # bonded mode return if (($mode>3)||($mode<0)); - + my $items = (2, 3, 4, 4)[$mode]; # items per entry my $name = ("bond", "angle", "dihedral", "improper")[$mode]; my $read = 0; @@ -691,7 +691,7 @@ } for (my $i=0; $i$last) { my $tmp = $first; $first = $last; $last = $tmp; } - if (($id2 = $hash{$hash_id = $first." ".$last}) eq "") + if (($id2 = $hash{$hash_id = $first." ".$last}) eq "") { $hash{$hash_id} = $id1; # add id to hash } @@ -757,26 +757,26 @@ } } - + sub AddMass { my $symbol = shift(@_); my $mass = shift(@_); - + return if ($symbols{$symbol} ne ""); $ids{++$max_id} = $symbol; $masses{$max_id} = $mass; $symbols{$symbol} = $max_id; } - + sub ReadTopology # read topology links { my $id = shift(@_); my $item = shift(@_); my $read = 0; my @tmp; - + open(TOPOLOGY, ") @@ -821,7 +821,7 @@ my @z = (-$L[2]/2, $L[2]/2); my $n = CRDGoto(atoms); my $extremes = !($L[0] && $L[1] && $L[2]); - + @Center = (0, 0, 0); return if (!$n); for (my $i=0; $i<$n; ++$i) @@ -858,7 +858,7 @@ my $s_OT = 1.7682; # CHARMM sigma [A] my $ahoh = (180-104.52)/360*PI(); my @p = ($loh*cos($ahoh), $loh*sin($ahoh), 0); - + printf("Info: creating fcc water\n") if ($info); $n_water = 4; # molecules/cell $nav = 6.022e23; # 1/mol @@ -867,11 +867,11 @@ @p_water = (0,0,0, @p, -$p[0],$p[1],0); $v_fcc = $n_water*$v_water; # cell volume $l_fcc = $v_fcc**(1/3); # cell length - @p_fcc = (0.00,0.00,0.00, 0.50,0.50,0.00, + @p_fcc = (0.00,0.00,0.00, 0.50,0.50,0.00, 0.50,0.00,0.50, 0.00,0.50,0.50); @n_fcc = (); for (my $i=0; $i0 ? int($x) : int($x)-1; } - + sub Periodic { my @p = splice(@_, 0, 3); - + return ( $p[0]-floor($p[0]/$L[0]+0.5)*$L[0], $p[1]-floor($p[1]/$L[1]+0.5)*$L[1], $p[2]-floor($p[2]/$L[2]+0.5)*$L[2]); } - - + + sub EraseWater { my $r = shift(@_)/2; @@ -918,7 +918,7 @@ $p[0]+$r,$p[1]+$r,$p[2]-$r, $p[0]+$r,$p[1]+$r,$p[2]+$r); my %list; my @n; - + my $d2 = ($r_water+$r)**2; my @l = ($L[0]/2, $L[1]/2, $L[2]/2); for (my $i=0; $i) - { - last if (!$n--); + while () + { + last if (!$n--); my @psf = split(" "); if ($res[1]!=$psf[2]) { ++$res[0]; $res[1] = $psf[2]; } printf(PSF_CTRL "%8d %4.4s %-4.4s %-4.4s %-4.4s %-4.4s ". @@ -1173,7 +1173,7 @@ my $n = PSFGoto(atoms); my $k = 0; my @res = (0, 0); - + CRDGoto(atoms); $net_charge = 0; printf(LAMMPS "Atoms%s\n\n",($add ? " # full" : "")) if ($n>0); @@ -1233,7 +1233,7 @@ { my @tmp = split(" "); printf(LAMMPS "%8d", ++$k); - for (my $j=0; $j<$n; ++$j) { + for (my $j=0; $j<$n; ++$j) { printf(LAMMPS " %16.12g", $j1)||!$water_dens); my $type = $mode ? CreateID(HT, OT, HT) : CreateID(HT, OT); my $id = $link{$type}; @@ -1266,7 +1266,7 @@ my $bit = 1; for (my $i=0; $i$project.in"); # use .in for temporary @@ -1448,7 +1448,7 @@ printf(LAMMPS "pair_style lj/charmm/coul/long 8 12\n"); printf(LAMMPS "pair_modify mix arithmetic\n"); printf(LAMMPS "kspace_style pppm 1e-6\n\n"); - + if ($cmap) { printf(LAMMPS "# Modify following line to point to the desired CMAP file\n"); printf(LAMMPS "fix cmap all cmap charmm$cmap.cmap\n"); @@ -1508,7 +1508,7 @@ # Molecular Dynamics Simulations, J. Comput. Chem. 25(2004): 1400-1415. # # ---------------------------------------------------------------------------- # - sub CharmmCmap + sub CharmmCmap { print "\nINITIATING CHARMM CMAP SUBROUTINE...\n\n"; @@ -1664,7 +1664,7 @@ # ... # # Criteria to be a PHI/PSI dihedral pair: - # 1. Atoms have to match with the mass/charge constellations as + # 1. Atoms have to match with the mass/charge constellations as # defined above. # 2. The atoms N--CA--C needs to be covalently bonded with each # other. @@ -1734,7 +1734,7 @@ my $N_PRO_counter = 0; my $C_flag = 0; - + for (my $i = 0; $i <= $natom_number; $i++) { my $cur_type = ${$atoms_matrix[$i]}[2]; my $cur_charge = ${$atoms_matrix[$i]}[3]; @@ -1777,8 +1777,8 @@ } # Quit if one of the atom types dosen't exist - if ( $C_counter == 0 or - ($CA_counter == 0 and $CA_GLY_counter == 0 and $CA_PRO_counter == 0) or + if ( $C_counter == 0 or + ($CA_counter == 0 and $CA_GLY_counter == 0 and $CA_PRO_counter == 0) or ($N_counter == 0 and $N_PRO_counter == 0) ) { if ($C_counter == 0) { print "\nCannot find the peptide backbone C atom type\n"; @@ -1814,21 +1814,21 @@ my $cur_atom4_type = ${atoms_matrix[${dihedrals_matrix[$i]}[5]-1]}[2]; next if (${dihedrals_matrix[$i]}[2] == ${dihedrals_matrix[$i-1]}[2] and - ${dihedrals_matrix[$i]}[3] == ${dihedrals_matrix[$i-1]}[3] and + ${dihedrals_matrix[$i]}[3] == ${dihedrals_matrix[$i-1]}[3] and ${dihedrals_matrix[$i]}[4] == ${dihedrals_matrix[$i-1]}[4] and ${dihedrals_matrix[$i]}[5] == ${dihedrals_matrix[$i-1]}[5]); # Determine PHI-dihedrals; If C-CA-N-C or C-N-CA-C, then save it in a list if ($cur_atom1_type == $C_type and $cur_atom4_type == $C_type) { - if ( ( ($cur_atom2_type == $CA_type or + if ( ( ($cur_atom2_type == $CA_type or $cur_atom2_type == $CA_GLY_type or - $cur_atom2_type == $CA_PRO_type) and + $cur_atom2_type == $CA_PRO_type) and ($cur_atom3_type == $N_type or - $cur_atom3_type == $N_PRO_type) ) or - ( ($cur_atom3_type == $CA_type or + $cur_atom3_type == $N_PRO_type) ) or + ( ($cur_atom3_type == $CA_type or $cur_atom3_type == $CA_GLY_type or - $cur_atom3_type == $CA_PRO_type) and - ($cur_atom2_type == $N_type or + $cur_atom3_type == $CA_PRO_type) and + ($cur_atom2_type == $N_type or $cur_atom2_type == $N_PRO_type) ) ) { push (@PHI_dihedrals,$cur_dihe_ID); $PHI_counter++; @@ -1837,17 +1837,17 @@ # Determin PSI-dihedrals; If N-CA-C-N or N-C-CA-N (N can be both normal N or N proline), # then save it in a list - if ( ($cur_atom1_type == $N_type and $cur_atom4_type == $N_type) or + if ( ($cur_atom1_type == $N_type and $cur_atom4_type == $N_type) or ($cur_atom4_type == $N_PRO_type and $cur_atom1_type == $N_PRO_type) or ($cur_atom1_type == $N_type and $cur_atom4_type == $N_PRO_type) or ($cur_atom4_type == $N_type and $cur_atom1_type == $N_PRO_type) ) { - if ( ( ($cur_atom2_type == $CA_type or + if ( ( ($cur_atom2_type == $CA_type or $cur_atom2_type == $CA_GLY_type or - $cur_atom2_type == $CA_PRO_type) and + $cur_atom2_type == $CA_PRO_type) and $cur_atom3_type == $C_type) or - ( ($cur_atom3_type == $CA_type or - $cur_atom3_type == $CA_GLY_type or - $cur_atom3_type == $CA_PRO_type) and + ( ($cur_atom3_type == $CA_type or + $cur_atom3_type == $CA_GLY_type or + $cur_atom3_type == $CA_PRO_type) and $cur_atom2_type == $C_type) ) { push (@PSI_dihedrals,$cur_dihe_ID); $PSI_counter++; @@ -1871,7 +1871,7 @@ # # The algorithm: # _____ - # | | + # | | # 1--2--3--4 PHI-dihedral # 4--3--2--1 # --C--N-CA--C--N-- Peptide backbone @@ -1887,7 +1887,7 @@ # if (2--3--4) = (4--3--2) # or # if (3--2--1) = (1--2--3) - # or + # or # if (3--2--1) = (4--3--2), # # then these 2 dihedrals are a PHI/PSI pair. If a pair is found, the @@ -1940,7 +1940,7 @@ if ($crossterm_CA_charge == $charge_CA) { $crossterm_type = 1; $crossterm_type1_flag = 1; } if ($crossterm_CA_charge == $charge_CA_GLY) { $crossterm_type = 5; $crossterm_type5_flag = 1; } - if ($crossterm_CA_charge == $charge_CA_PRO) { + if ($crossterm_CA_charge == $charge_CA_PRO) { $crossterm_type = 3; $crossterm_type3_flag = 1; # Checking the last crossterm, re-assign the last crossterm type if needed if ($crossterm_counter-1 >= 0 and $PHI_PSI_matrix[$crossterm_counter-1][0] == 1) { @@ -1998,7 +1998,7 @@ print "acid abbreviation for that amino acid, and <#2> is the molecule number\n"; print "of the terminal amino acid residue.\n\n"; print "For example, if the last atom of the last amino acid in the peptide\n"; - print "sequence is listed in the pdb file as:\n\n"; + print "sequence is listed in the pdb file as:\n\n"; print " 'ATOM 853 O GLU P 56 12.089 -1.695 -6.543 1.00 1.03 PROA'\n\n"; print "you would insert the following line after it:\n\n"; print " 'TER 854 GLU 56'\n\n"; @@ -2016,7 +2016,7 @@ # Print out PHI/PSI diheral pair list my $pair_counter = 0; # Don't presently use $ncrosstermtypes but have this available if wish to print it out - my $ncrosstermtypes = $crossterm_type1_flag + $crossterm_type2_flag + $crossterm_type3_flag + + my $ncrosstermtypes = $crossterm_type1_flag + $crossterm_type2_flag + $crossterm_type3_flag + $crossterm_type4_flag + $crossterm_type5_flag + $crossterm_type6_flag; print "\nWriting \"$project.data\" with section \"CMAP crossterms\" added at the end.\n"; diff --git a/tools/ch2lmp/lammps2pdb.pl b/tools/ch2lmp/lammps2pdb.pl index 8bf636ea4f..f1180023a5 100755 --- a/tools/ch2lmp/lammps2pdb.pl +++ b/tools/ch2lmp/lammps2pdb.pl @@ -8,9 +8,9 @@ # conjunction with vmd # # Notes: Copyright by author for Sandia National Laboratories -# 20040903 Conception date of v1.0: rudimentary script for collagen +# 20040903 Conception date of v1.0: rudimentary script for collagen # project. -# 20050423 Conception date of v2.0: +# 20050423 Conception date of v2.0: # - changed data access through indexing data directly on disk; # - added all command line options # 20050425 Corrected focussing to use a target molecule's moment of @@ -22,21 +22,21 @@ # the data stream. # subroutines - + sub test { my $name = shift(@_); - + printf("Error: file %s not found\n", $name) if (!scalar(stat($name))); return !scalar(stat($name)); } - + sub initialize { my $k = 0; - my @options = ("-help", "-nstart", "-dn", "-cut", "-repair", - "-units", "-quiet", "-nodetect", "-data", "-pbc", + my @options = ("-help", "-nstart", "-dn", "-cut", "-repair", + "-units", "-quiet", "-nodetect", "-data", "-pbc", "-focus", "-center", "-exclude"); my @remarks = ("display this message", "starting frame [-1]", @@ -61,7 +61,7 @@ "* Expected files in current directory: project.data, project.dump", "* Generated output files: project_trj.psf, project_trj.pdb", "* Uses project_ctrl.psf if available"); - + $program = "lammps2pdb"; $version = "2.2.5"; @@ -96,7 +96,7 @@ $info = $switch ? 0 : 1 if (!$k--); $detect = $switch ? 0 : 1 if (!$k--); $data_name = $arg[1] if (!$k--); - if (!$k--) { + if (!$k--) { if ($switch) { $pbc{ALL} = 1; } else { foreach (split(",",$arg[1])) { $pbc{uc($_)} = 1; }}} if (!$k--) { foreach (split(",",$arg[1])) { $focus{uc($_)} = uc($_);}} @@ -120,7 +120,7 @@ printf("\n"); exit(-1); } - printf("%s v%s (c)%s\n\n", $program, $version, $year) if ($info); + printf("%s v%s (c)%s\n\n", $program, $version, $year) if ($info); $data_name = $project.".data" if ($data_name eq ""); $traject_name = $project.".dump"; @@ -133,7 +133,7 @@ my $flag = test($data_name); printf("Conversion aborted\n\n") if ($flag); exit(-1) if ($flag); - + # data input create_atom_ids(); @@ -145,9 +145,9 @@ open(PSF, ">".$psf_name) if (!$psf_ctrl); open(PDB, ">".$pdb_name); - + # align center with focus - + %center = %focus if (scalar(%focus)); } @@ -166,24 +166,24 @@ return @_; } - + sub V_Add # V_Add(@a, @b) = @a + @b; { return (@_[0]+@_[3], @_[1]+@_[4], @_[2]+@_[5]); } - - + + sub V_Subtr # V_Subtr(@a, @b) = @a - @b; { return (@_[0]-@_[3], @_[1]-@_[4], @_[2]-@_[5]); } - - + + sub V_Dot # V_Dot(@a, @b) = @a . @b; { return (@_[0]*@_[3]+@_[1]*@_[4]+@_[2]*@_[5]); } - + sub V_Cross # V_Cross(@a, @b) = @a x @b; { @@ -191,7 +191,7 @@ @_[0]*@_[4]-@_[1]*@_[3]); } - + sub V_Mult # V_Mult($f, @a) = $f * @a; { return (@_[0]*@_[1], @_[0]*@_[2], @_[0]*@_[3]); @@ -203,11 +203,11 @@ return V_Mult(1/sqrt(V_Dot(@_[0,1,2],@_[0,1,2])), @_[0,1,2]); } - + sub pbc # periodic -0.5*L <= x < 0.5*L { my $x = @_[0]/@_[1]+0.5; - + return @_[0]-@_[1]*($x<0 ? int($x)-1 : int($x)); } @@ -216,8 +216,8 @@ { return (pbc(@_[0], @_[3]), pbc(@_[1], @_[4]), pbc(@_[2], @_[5])); } - - + + sub V_Cmp # V_Cmp(abs(@a), abs(@b)) { return -1 if (abs($_[0])abs($_[5])); return 0; } - - + + sub V_Sort # sort on descending absolute { # value my @v = @_; - + for (my $i=0; $i 0 @A = ($f2<0 ? abs(2*$f2)**(1/3) : 0, 0); } - else { + else { if (abs($f2)<1e-14) { # limit f2 -> 0 my $f = sqrt(abs($f1))/$c3; @A = $f1<0 ? (-$f*$z1[1], 0.5*$f) : ($f, 0); @B = $f1<0 ? (-$A[0], $A[1]) : ($f, 0); } else { - @B = C_Pow(C_Add(($f2, 0), + @B = C_Pow(C_Add(($f2, 0), C_Pow(($f1*$f1*$f1+$f2*$f2, 0), 1/2)), 1/3); @A = C_Div(($f1/$c3, 0), @B); @B = ($B[0]/$c3, $B[1]/$c3); } } @@ -503,21 +503,21 @@ my $input = shift; my $dlines = shift; my $read; - + while ($dlines--) { $read = <$input>; } return $read; } - + sub rewind_read { my $input = shift; - + seek($input, 0, SEEK_SET); } -# create crossreference tables +# create crossreference tables sub create_psf_index # make an index of id { # locations @@ -538,13 +538,13 @@ } } - + sub psf_goto # goto $ident in { create_psf_index() if (!scalar(%PSFIndex)); my $id = shift(@_); my @n = split(" ", $PSFIndex{$id}); - + @PSFBuffer = (); if (!scalar(@n)) { @@ -571,7 +571,7 @@ my $id; my %hash; my %size; - + foreach ((masses,atoms,bonds,angles,dihedrals,impropers)) { $hash{$_}=$_; } open(DATA, "<".$data_name); for (my $i=0; $i<2; ++$i) { my $tmp = ; } # skip first two lines @@ -604,7 +604,7 @@ } } - + sub goto_data { create_data_index() if (!scalar(%DATAIndex)); @@ -623,7 +623,7 @@ # create atom and residue identifiers - + sub create_names { return if (scalar(@names)); @@ -636,7 +636,7 @@ my @letter = ("X", "Y", "Z", "P", "Q", "R", "S", "A", "B", "C"); my $k = 0; my %atom; - + $names[0] = ""; foreach (@mass) { $atom{$_} = shift(@name); } for (my $i=1; $i<=$n; ++$i) @@ -658,7 +658,7 @@ my @data = @_; my $p = $data[1]." ".$data[2]; my $k; - + for ($k=0; ($k))), 0,length($subject)) ne $subject)) {} } - - + + sub read_traject { my @box; my @l; - + advance_traject("timestep"); my $timestep = (split(" ", ))[0]; advance_traject("number of atoms"); my $n = (split(" ", ))[0]; advance_traject("box bounds"); for (my $i=0; $i<3; ++$i) - { + { my @data = split(" ", ); $box[$i] = $data[0]; # box edge $l[$i] = $data[1]-$data[0]; # box length @@ -806,9 +806,9 @@ return ($timestep, $n, @box, @l); } - + # pdb format - + sub eigen_vector # eigen_vector(@A, $l) { my @A = splice(@_,0,9); @@ -826,8 +826,8 @@ return (0,0,1) if ($A[8]==1); return (0,0,0); } - - + + sub pdb_inertia { my @s = ( @@ -852,13 +852,13 @@ M_Transpose(M_RotateNormal(MV_Dot(@A,@b)))); return M_Dot(@B, @A); } - - + + sub pdb_focus # using moment of inertia { my @R = pdb_inertia(@_); - + printf("Info: focussing\n") if ($info); foreach (@position) { @@ -867,7 +867,7 @@ } } - + sub pdb_center { my @c = splice(@_, 0, 3); @@ -881,7 +881,7 @@ } } - + sub pdb_pbc { printf("Info: applying periodicity\n") if ($info); @@ -909,7 +909,7 @@ } } - + sub pdb_positions { my @m = (0,0,0,0,0,0,0,0,0); @@ -918,7 +918,7 @@ my $mass; my @p; my $d; - + foreach (@traject) { my @arg = split(" "); @@ -974,7 +974,7 @@ printf(PDB "%-11.11s", "P 1"); printf(PDB "%3.3s\n", "1"); } - + sub pdb_atoms { @@ -998,15 +998,15 @@ foreach (@p) { $_ = 0 if (abs($_)<1e-4); } printf(PDB "ATOM "); # pdb command printf(PDB "%6.6s ", ++$n); # atom number - printf(PDB "%-3.3s ", + printf(PDB "%-3.3s ", $psf_ctrl ? $psf[4] : $names[$p[4]]); # atom name - printf(PDB "%-3.3s ", + printf(PDB "%-3.3s ", $psf_ctrl ? $psf[3] : $residue[$nres]); # residue name printf(PDB "%5.5s ", $nres); # residue number printf(PDB "%3.3s ", ""); # empty placeholder printf(PDB "%7.7s %7.7s %7.7s ", $p[0], $p[1], $p[2]); # positions - printf(PDB "%5.5s %5.5s %4.4s ", + printf(PDB "%5.5s %5.5s %4.4s ", "1.00", "0.00", ""); # trailing variables printf(PDB "%-4.4s\n", $psf_ctrl ? $psf[1] : $cluster[$nres]); # cluster name @@ -1015,7 +1015,7 @@ }; printf(PDB "END\n"); } - + sub pdb_timestep { @@ -1045,7 +1045,7 @@ my $l = 0; my $k = 0; my @extra; - + for (my $i=0; $i<$n; ++$i) { my @arg = split(" ", ); @@ -1073,13 +1073,13 @@ } printf(PSF "\n"); } - + sub psf_bonds { my $npairs = 0; - - delete_exclude() if (scalar(%exclude)>0); + + delete_exclude() if (scalar(%exclude)>0); delete_crossovers() if ($cut); printf(PSF "%8.8s !NBOND\n", scalar(@bonds)); foreach (@bonds) @@ -1098,7 +1098,7 @@ initialize(); # create .pdb file - + $ncurrent = -1; while ($traject_flag&&!eof(TRAJECT)) { @@ -1134,4 +1134,4 @@ close(TRAJECT); close(DATA); - +