1 #!/bin/env python 2 # 3 # File: PyMOLCalculatePhiPsiAngles.py 4 # Author: Manish Sud <msud@san.rr.com> 5 # 6 # Copyright (C) 2026 Manish Sud. All rights reserved. 7 # 8 # The functionality available in this script is implemented using PyMOL, a 9 # molecular visualization system on an open source foundation originally 10 # developed by Warren DeLano. 11 # 12 # This file is part of MayaChemTools. 13 # 14 # MayaChemTools is free software; you can redistribute it and/or modify it under 15 # the terms of the GNU Lesser General Public License as published by the Free 16 # Software Foundation; either version 3 of the License, or (at your option) any 17 # later version. 18 # 19 # MayaChemTools is distributed in the hope that it will be useful, but without 20 # any warranty; without even the implied warranty of merchantability of fitness 21 # for a particular purpose. See the GNU Lesser General Public License for more 22 # details. 23 # 24 # You should have received a copy of the GNU Lesser General Public License 25 # along with MayaChemTools; if not, see <http://www.gnu.org/licenses/> or 26 # write to the Free Software Foundation Inc., 59 Temple Place, Suite 330, 27 # Boston, MA, 02111-1307, USA. 28 # 29 30 from __future__ import print_function 31 32 import os 33 import sys 34 import time 35 import re 36 37 # PyMOL imports... 38 try: 39 import pymol 40 41 # Finish launching PyMOL in a command line mode for batch processing (-c) 42 # along with the following options: disable loading of pymolrc and plugins (-k); 43 # suppress start up messages (-q) 44 pymol.finish_launching(["pymol", "-ckq"]) 45 except ImportError as ErrMsg: 46 sys.stderr.write("\nFailed to import PyMOL module/package: %s\n" % ErrMsg) 47 sys.stderr.write("Check/update your PyMOL environment and try again.\n\n") 48 sys.exit(1) 49 50 # MayaChemTools imports... 51 sys.path.insert(0, os.path.join(os.path.dirname(sys.argv[0]), "..", "lib", "Python")) 52 try: 53 from docopt import docopt 54 import MiscUtil 55 import PyMOLUtil 56 except ImportError as ErrMsg: 57 sys.stderr.write("\nFailed to import MayaChemTools module/package: %s\n" % ErrMsg) 58 sys.stderr.write("Check/update your MayaChemTools environment and try again.\n\n") 59 sys.exit(1) 60 61 ScriptName = os.path.basename(sys.argv[0]) 62 Options = {} 63 OptionsInfo = {} 64 65 66 def main(): 67 """Start execution of the script.""" 68 69 MiscUtil.PrintInfo( 70 "\n%s (PyMOL v%s; MayaChemTools v%s; %s): Starting...\n" 71 % (ScriptName, pymol.cmd.get_version()[0], MiscUtil.GetMayaChemToolsVersion(), time.asctime()) 72 ) 73 74 (WallClockTime, ProcessorTime) = MiscUtil.GetWallClockAndProcessorTime() 75 76 # Retrieve command line arguments and options... 77 RetrieveOptions() 78 79 # Process and validate command line arguments and options... 80 ProcessOptions() 81 82 # Perform actions required by the script... 83 CalculatePhiPsiAngles() 84 85 MiscUtil.PrintInfo("\n%s: Done...\n" % ScriptName) 86 MiscUtil.PrintInfo("Total time: %s" % MiscUtil.GetFormattedElapsedTime(WallClockTime, ProcessorTime)) 87 88 89 def CalculatePhiPsiAngles(): 90 """Calculate phi and psi angles for macromolecules containing amino acids.""" 91 92 SetupOutputFiles() 93 WriteColumnLabels() 94 95 Infile = OptionsInfo["Infile"] 96 MolName = OptionsInfo["InfileRoot"] 97 98 MiscUtil.PrintInfo("\nCalculating phi and psi torsion angles for input file %s..." % Infile) 99 100 # Load infile 101 pymol.cmd.load(Infile, MolName) 102 103 OutDelim = OptionsInfo["OutDelim"] 104 Precision = OptionsInfo["Precision"] 105 106 # Go over specified chain IDs.. 107 for ChainID in OptionsInfo["SpecifiedChainsAndLigandsInfo"]["ChainIDs"]: 108 # Write out information for combined file... 109 PhiPsiInfo = PyMOLUtil.GetPhiPsiResiduesInfo(MolName, ChainID, Categorize=True) 110 OptionsInfo["OutfileResCount"] += len(PhiPsiInfo["ResNums"]) 111 WritePhiPsiInfo(OptionsInfo["OutFH"], MolName, ChainID, PhiPsiInfo, OutDelim, Precision) 112 113 # Write out information for category fies... 114 if OptionsInfo["MultipleOutFiles"]: 115 PhiPsiInfoList = [] 116 GeneralPhiPsiInfo, GlycinePhiPsiInfo, ProlinePhiPsiInfo, PreProlinePhiPsiInfo = ( 117 PyMOLUtil.GetPhiPsiCategoriesResiduesInfo(MolName, ChainID) 118 ) 119 PhiPsiInfoList.extend([GeneralPhiPsiInfo, GlycinePhiPsiInfo, ProlinePhiPsiInfo, PreProlinePhiPsiInfo]) 120 121 for Index, Category in enumerate(OptionsInfo["Categories"]): 122 OptionsInfo["CategoriesResCount"][Category] += len(PhiPsiInfoList[Index]["ResNums"]) 123 WritePhiPsiInfo( 124 OptionsInfo["CategoriesOutFHs"][Category], 125 MolName, 126 ChainID, 127 PhiPsiInfoList[Index], 128 OutDelim, 129 Precision, 130 ) 131 132 # Delete MolName object 133 pymol.cmd.delete(MolName) 134 135 # Close all files... 136 CloseOutputFiles() 137 138 # List number of phi and psi angles in output files... 139 MiscUtil.PrintInfo( 140 "\nNumber of phi and psi angles in output file %s: %d" 141 % (OptionsInfo["Outfile"], OptionsInfo["OutfileResCount"]) 142 ) 143 if OptionsInfo["MultipleOutFiles"]: 144 MiscUtil.PrintInfo("") 145 for Index, Category in enumerate(OptionsInfo["Categories"]): 146 MiscUtil.PrintInfo( 147 "Number of phi and psi angles in output file %s: %d" 148 % (OptionsInfo["CategoriesOutfiles"][Category], OptionsInfo["CategoriesResCount"][Category]) 149 ) 150 151 152 def WritePhiPsiInfo(OutFH, MolName, ChainID, PhiPsiResiduesInfo, OutDelim, Precision): 153 """Write out phi and psi information.""" 154 155 for ResNum in PhiPsiResiduesInfo["ResNums"]: 156 ResName = PhiPsiResiduesInfo["ResName"][ResNum] 157 Phi = "%.*f" % (Precision, PhiPsiResiduesInfo["Phi"][ResNum]) 158 Psi = "%.*f" % (Precision, PhiPsiResiduesInfo["Psi"][ResNum]) 159 Category = PhiPsiResiduesInfo["Category"][ResNum] 160 161 LineWords = [] 162 if OptionsInfo["OutChainID"]: 163 LineWords.append(ChainID) 164 165 LineWords.extend([ResNum, ResName, Phi, Psi]) 166 167 if OptionsInfo["OutCategory"]: 168 LineWords.append(Category) 169 170 Line = OutDelim.join(LineWords) 171 172 OutFH.write("%s\n" % Line) 173 174 175 def WriteColumnLabels(): 176 """Write out column labels.""" 177 178 OutDelim = OptionsInfo["OutDelim"] 179 180 ColLabels = [] 181 if OptionsInfo["OutChainID"]: 182 ColLabels.append("ChainID") 183 184 ColLabels.extend(["ResNum", "ResName", "Phi", "Psi"]) 185 186 if OptionsInfo["OutCategory"]: 187 ColLabels.append("Category") 188 189 Line = OutDelim.join(ColLabels) 190 191 # Write column labels for combined file... 192 OutFH = OptionsInfo["OutFH"] 193 OutFH.write("%s\n" % Line) 194 195 if not OptionsInfo["MultipleOutFiles"]: 196 return 197 198 # Write columns labels for category files... 199 for Category in OptionsInfo["Categories"]: 200 CategoryOutFH = OptionsInfo["CategoriesOutFHs"][Category] 201 CategoryOutFH.write("%s\n" % Line) 202 203 204 def SetupOutputFiles(): 205 """Open output files.""" 206 207 if OptionsInfo["MultipleOutFiles"]: 208 MiscUtil.PrintInfo("\nGenerating output files: %s" % (", ".join(OptionsInfo["OutfilesList"]))) 209 else: 210 MiscUtil.PrintInfo("\nGenerating output file %s..." % (OptionsInfo["Outfile"])) 211 212 # Open combined output file... 213 Outfile = OptionsInfo["Outfile"] 214 OutFH = open(Outfile, "w") 215 if OutFH is None: 216 MiscUtil.PrintError("Couldn't open output file: %s.\n" % (Outfile)) 217 OptionsInfo["OutFH"] = OutFH 218 OptionsInfo["OutfileResCount"] = 0 219 220 if not OptionsInfo["MultipleOutFiles"]: 221 return 222 223 # Open output files for different categories... 224 OptionsInfo["CategoriesOutFHs"] = {} 225 OptionsInfo["CategoriesResCount"] = {} 226 for Category in OptionsInfo["Categories"]: 227 CategoryOutfile = OptionsInfo["CategoriesOutfiles"][Category] 228 CategoryOutFH = open(CategoryOutfile, "w") 229 if CategoryOutfile is None: 230 MiscUtil.PrintError("Couldn't open output file: %s.\n" % (CategoryOutfile)) 231 232 OptionsInfo["CategoriesOutFHs"][Category] = CategoryOutFH 233 OptionsInfo["CategoriesResCount"][Category] = 0 234 235 236 def CloseOutputFiles(): 237 """Close output files.""" 238 239 OptionsInfo["OutFH"].close() 240 241 if not OptionsInfo["MultipleOutFiles"]: 242 return 243 244 for Category in OptionsInfo["Categories"]: 245 CategoryOutFH = OptionsInfo["CategoriesOutFHs"][Category] 246 CategoryOutFH.close() 247 248 249 def RetrieveInfileInfo(): 250 """Retrieve information for input file.""" 251 252 Infile = OptionsInfo["Infile"] 253 InfileRoot = OptionsInfo["InfileRoot"] 254 255 ChainsAndLigandsInfo = PyMOLUtil.GetChainsAndLigandsInfo(Infile, InfileRoot) 256 OptionsInfo["ChainsAndLigandsInfo"] = ChainsAndLigandsInfo 257 258 259 def ProcessChainIDs(): 260 """Process specified chain IDs for infile.""" 261 262 MiscUtil.PrintInfo("\nProcessing specified chain IDs for input file %s..." % OptionsInfo["Infile"]) 263 ChainsAndLigandsInfo = OptionsInfo["ChainsAndLigandsInfo"] 264 SpecifiedChainsAndLigandsInfo = PyMOLUtil.ProcessChainsAndLigandsOptionsInfo( 265 ChainsAndLigandsInfo, "-c, --chainIDs", OptionsInfo["ChainIDs"], None, None 266 ) 267 268 OptionsInfo["SpecifiedChainsAndLigandsInfo"] = SpecifiedChainsAndLigandsInfo 269 270 MiscUtil.PrintInfo("Specified chain IDs: %s" % (", ".join(SpecifiedChainsAndLigandsInfo["ChainIDs"]))) 271 272 273 def SetupCategoryOutfiles(): 274 """Setup output file names for different categories of phi and psi angles.""" 275 276 # Initialize... 277 OptionsInfo["OutfilesList"] = [] 278 OptionsInfo["OutfilesList"].append(OptionsInfo["Outfile"]) 279 280 OptionsInfo["Categories"] = ["General", "Glycine", "Proline", "PreProline"] 281 OptionsInfo["CategoriesOutfiles"] = {} 282 for Category in OptionsInfo["Categories"]: 283 OptionsInfo["CategoriesOutfiles"][Category] = None 284 285 if not OptionsInfo["MultipleOutFiles"]: 286 return 287 288 FileDir, FileName, FileExt = MiscUtil.ParseFileName(OptionsInfo["Outfile"]) 289 OutfileRoot = FileName 290 OutfileExt = FileExt 291 292 for Category in OptionsInfo["Categories"]: 293 CategoryOutfile = "%s_%s.%s" % (OutfileRoot, Category, OutfileExt) 294 if os.path.exists(CategoryOutfile): 295 if not OptionsInfo["Overwrite"]: 296 MiscUtil.PrintError( 297 'The category output file, %s, already exist. Use option "--ov" or "--overwrite" and try again.\n' 298 % (CategoryOutfile) 299 ) 300 301 OptionsInfo["CategoriesOutfiles"][Category] = CategoryOutfile 302 OptionsInfo["OutfilesList"].append(CategoryOutfile) 303 304 305 def ProcessOptions(): 306 """Process and validate command line arguments and options.""" 307 308 MiscUtil.PrintInfo("Processing options...") 309 310 # Validate options... 311 ValidateOptions() 312 313 OptionsInfo["Infile"] = Options["--infile"] 314 FileDir, FileName, FileExt = MiscUtil.ParseFileName(OptionsInfo["Infile"]) 315 OptionsInfo["InfileRoot"] = FileName 316 317 OptionsInfo["Outfile"] = Options["--outfile"] 318 FileDir, FileName, FileExt = MiscUtil.ParseFileName(OptionsInfo["Outfile"]) 319 OptionsInfo["OutfileRoot"] = FileName 320 321 OptionsInfo["Overwrite"] = Options["--overwrite"] 322 323 OptionsInfo["OutDelim"] = " " 324 if MiscUtil.CheckFileExt(OptionsInfo["Outfile"], "csv"): 325 OptionsInfo["OutDelim"] = "," 326 elif MiscUtil.CheckFileExt(OptionsInfo["Outfile"], "tsv txt"): 327 OptionsInfo["OutDelim"] = "\t" 328 else: 329 MiscUtil.PrintError( 330 'The file name specified , %s, for option "--outfile" is not valid. Supported file formats: csv tsv txt\n' 331 % (OptionsInfo["Outfile"]) 332 ) 333 334 OptionsInfo["OutMode"] = Options["--outMode"] 335 OptionsInfo["MultipleOutFiles"] = True if re.match("^MultipleFiles$", OptionsInfo["OutMode"], re.I) else False 336 337 OptionsInfo["OutChainID"] = True if re.match("^Yes$", Options["--outChainID"], re.I) else False 338 OptionsInfo["OutCategory"] = True if re.match("^Yes$", Options["--outCategory"], re.I) else False 339 340 OptionsInfo["Overwrite"] = Options["--overwrite"] 341 OptionsInfo["Precision"] = int(Options["--precision"]) 342 343 RetrieveInfileInfo() 344 345 OptionsInfo["ChainIDs"] = Options["--chainIDs"] 346 ProcessChainIDs() 347 348 SetupCategoryOutfiles() 349 350 351 def RetrieveOptions(): 352 """Retrieve command line arguments and options.""" 353 354 # Get options... 355 global Options 356 Options = docopt(_docoptUsage_) 357 358 # Set current working directory to the specified directory... 359 WorkingDir = Options["--workingdir"] 360 if WorkingDir: 361 os.chdir(WorkingDir) 362 363 # Handle examples option... 364 if "--examples" in Options and Options["--examples"]: 365 MiscUtil.PrintInfo(MiscUtil.GetExamplesTextFromDocOptText(_docoptUsage_)) 366 sys.exit(0) 367 368 369 def ValidateOptions(): 370 """Validate option values.""" 371 372 MiscUtil.ValidateOptionFilePath("-i, --infile", Options["--infile"]) 373 MiscUtil.ValidateOptionFileExt("-i, --infile", Options["--infile"], "pdb cif") 374 375 MiscUtil.ValidateOptionFileExt("-o, --outfile", Options["--outfile"], "csv tsv txt") 376 MiscUtil.ValidateOptionsOutputFileOverwrite( 377 "-o, --outfile", Options["--outfile"], "--overwrite", Options["--overwrite"] 378 ) 379 380 MiscUtil.ValidateOptionTextValue("--outMode", Options["--outMode"], "SingleFile MultipleFiles") 381 MiscUtil.ValidateOptionTextValue("--outChainID", Options["--outChainID"], "yes no") 382 MiscUtil.ValidateOptionTextValue("--outCategory", Options["--outCategory"], "yes no") 383 MiscUtil.ValidateOptionIntegerValue("-p, --precision", Options["--precision"], {">": 0}) 384 385 386 # Setup a usage string for docopt... 387 _docoptUsage_ = """ 388 PyMOLCalculatePhiPsiAngles.py - Calculate phi and psi torsion angles 389 390 Usage: 391 PyMOLCalculatePhiPsiAngles.py [--chainIDs <First, All or ID1,ID2...>] 392 [--outMode <SingleFile or MultipleFies>] [--outChainID <yes or no>] 393 [--outCategory <yes or no>] [--overwrite] [--precision <number>] 394 [-w <dir>] -i <infile> -o <outfile> 395 PyMOLCalculatePhiPsiAngles.py -h | --help | -e | --examples 396 397 Description: 398 Calculate phi and psi torsion angels for amino acid residues present 399 in macromolecules. 400 401 The phi and psi angles are categorized into the following groups 402 corresponding to four types of Ramachandran plots: 403 404 General: All residues except glycine, proline, or pre-proline 405 Glycine: Only glycine residues 406 Proline: Only proline residues 407 Pre-Proline: Only residues before proline not including glycine or 408 proline 409 410 The supported input file format are: PDB (.pdb), mmCIF (.cif) 411 412 The supported output file formats are: CSV (.csv), TSV (.tsv, .txt) 413 414 Options: 415 -c, --chainIDs <First, All or ID1,ID2...> [default: All] 416 List of chain IDs to use for calculating phi and psi angles for residues 417 in chains. Possible values: First, All, or a comma delimited list of chain 418 IDs. The default is to use all chain IDs in input file. 419 -e, --examples 420 Print examples. 421 -h, --help 422 Print this help message. 423 -i, --infile <infile> 424 Input file name. 425 -o, --outfile <outfile> 426 Output file name for writing out calculated values. Supported text file 427 extensions: csv, tsv or txt. 428 429 In addition to the specified outfile containing phi and psi angles for all 430 residues, a set of additional output files is generated for 'MultipleFiles' 431 value of '--outMode' option. The names of these output files are 432 automatically generated from the the name of the specified output 433 file as shown below: 434 435 General: <OutfileRoot>_General.<OutfileExt> 436 Glycine: <OutfileRoot>_Glycine.<OutfileExt> 437 Proline: <OutfileRoot>_Proline.<OutfileExt> 438 Pre-Proline: <OutfileRoot>_PreProline.<OutfileExt> 439 440 --outMode <SingleFile or MultipleFiles> [default: SingleFile] 441 A single output file containing phi and psi angles for all residues or 442 multiple output files corresponding to different categories of angles. 443 444 The phi and psi angles are categorized into the following groups 445 corresponding to four types of Ramachandran plots: 446 447 General: All residues except glycine, proline, or pre-proline 448 Glycine: Only glycine residues 449 Proline: Only proline residues 450 Pre-Proline: Only residues before proline not including glycine or 451 proline 452 453 The output files contain the following information: 454 455 ChainID ResNum ResName Phi Psi Category 456 457 --outChainID <yes or no> [default: yes] 458 Write chain IDs to output file. 459 --outCategory <yes or no> [default: yes] 460 Write phi and psi category to output file. 461 --overwrite 462 Overwrite existing files. 463 -p, --precision <number> [default: 2] 464 Floating point precision for writing the calculated phi and psi angles. 465 -w, --workingdir <dir> 466 Location of working directory which defaults to the current directory. 467 468 Examples: 469 To calculate phi and psi angles for all residues across all chains in input 470 file and write out a single CSV file containing calculated values along with 471 chain IDs, residue names and numbers, and category of angles corresponding 472 to Ramachandran plots, type: 473 474 % PyMOLCalculatePhiPsiAngles.py -i Sample3.pdb -o Sample3Out.csv 475 476 To calculate phi and psi angles for all residues across all chains in input 477 file and write out a multiple CSV files corresponding to categories of angles 478 for Ramachandran plots along with other relevant information, type: 479 480 % PyMOLCalculatePhiPsiAngles.py --outMode MultipleFiles -i Sample3.pdb 481 -o Sample3Out.csv 482 483 To calculate phi and psi angles for all residues in a specific chain in input 484 file and write out a single TSV file containing calculated values along with 485 other relevant information, type: 486 487 % PyMOLCalculatePhiPsiAngles.py -c E -i Sample3.pdb -o Sample3Out.csv 488 489 To calculate phi and psi angles for all residues in a specific chain in input 490 file and write out a multiple TSV files containing calculated values at a specific 491 precision along with other relevant information, type: 492 493 % PyMOLCalculatePhiPsiAngles.py --outMode MultipleFiles --chainIDs I 494 -i Sample3.pdb -o Sample3Out.csv 495 496 Author: 497 Manish Sud(msud@san.rr.com) 498 499 See also: 500 DownloadPDBFiles.pl, PyMOLCalculateRMSD.py, PyMOLCalculateProperties.py, 501 PyMOLGenerateRamachandranPlots.py 502 503 Copyright: 504 Copyright (C) 2026 Manish Sud. All rights reserved. 505 506 The functionality available in this script is implemented using PyMOL, a 507 molecular visualization system on an open source foundation originally 508 developed by Warren DeLano. 509 510 This file is part of MayaChemTools. 511 512 MayaChemTools is free software; you can redistribute it and/or modify it under 513 the terms of the GNU Lesser General Public License as published by the Free 514 Software Foundation; either version 3 of the License, or (at your option) any 515 later version. 516 517 """ 518 519 if __name__ == "__main__": 520 main()
