#!/usr/bin/env python3 # (C) 2023 Cadence Design Systems, Inc. (Cadence) # All rights reserved. # TERMS FOR USE OF SAMPLE CODE The software below ("Sample Code") is # provided to current licensees or subscribers of Cadence products or # SaaS offerings (each a "Customer"). # Customer is hereby permitted to use, copy, and modify the Sample Code, # subject to these terms. Cadence claims no rights to Customer's # modifications. Modification of Sample Code is at Customer's sole and # exclusive risk. Sample Code may require Customer to have a then # current license or subscription to the applicable Cadence offering. # THE SAMPLE CODE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, # EXPRESS OR IMPLIED. OPENEYE DISCLAIMS ALL WARRANTIES, INCLUDING, BUT # NOT LIMITED TO, WARRANTIES OF MERCHANTABILITY, FITNESS FOR A # PARTICULAR PURPOSE AND NONINFRINGEMENT. In no event shall Cadence be # liable for any damages or liability in connection with the Sample Code # or its use. """Depict the B-factor of a ligand and its environment.""" import argparse import io import os import sys from pathlib import Path from openeye import oechem, oedepict, oegrapheme from PIL import Image from rich_argparse import HelpPreviewAction, RichHelpFormatter __SCRIPT_NAME__ = Path(__file__).absolute().stem __SCRIPT_DESC__ = "Depict B-factor of an active site." __SCRIPT_TOOLKITS__ = ["oechem", "oedepict", "oegrapheme"] __SCRIPT_CATEGORIES__ = ["visualization"] def parse_options() -> argparse.Namespace: """Set up command line options.""" parser = argparse.ArgumentParser( add_help=True, formatter_class=RichHelpFormatter, description="[yellow]" + __SCRIPT_DESC__ + "[/yellow]", ) parser.add_argument("--help-image", action=HelpPreviewAction) # input options input_group = parser.add_argument_group("Input ligand-protein complex") exclusive_input_group = input_group.add_mutually_exclusive_group(required=True) exclusive_input_group.add_argument( "--complex", type=str, required=False, metavar="PDB-FILE", help="input PDB file of the ligand-protein complex", ) exclusive_input_group.add_argument( "--design-unit", "--du", type=str, metavar="DU-FILE", help="input design unit file", ) image_group = parser.add_argument_group("Image options") image_group.add_argument( "--image", type=str, required=False, metavar="IMAGE-FILE", help="output image file (PNG, SVG) (required: %(required)s) -- if no output is provided the image will be displayed on the screen", ) image_group.add_argument( "--width", type=int, default=600, help="width of output image (default: %(default)s)", ) image_group.add_argument( "--height", type=int, default=400, help="height of output image (default: %(default)s)", ) image_group.add_argument( "--interactive-legend", default=False, action="store_true", help="visualize legend on mouse hover (SVG-only feature) (default: %(default)s)", ) return parser.parse_args() def main() -> int: """Depict B-factor.""" args = parse_options() _check_image_file(args) if args.complex: protein, ligand = get_protein_and_ligand_from_pdb(args.complex) elif args.design_unit: protein, ligand = get_protein_and_ligand_from_design_unit(args.design_unit) else: oechem.OEThrow.Fatal("Invalid input option!") # depict active site with b-bfactor image = oedepict.OEImage(args.width, args.height) opts = oegrapheme.OE2DActiveSiteDisplayOptions(args.width, args.height) opts.SetRenderInteractiveLegend(args.interactive_legend) depict_bfactor_map(image, protein, ligand, opts) if args.image and Path(args.image).suffix[1:].lower() == "svg": icon_scale = 0.5 oedepict.OEAddInteractiveIcon( image, oedepict.OEIconLocation_TopRight, icon_scale ) oedepict.OEDrawCurvedBorder(image, oedepict.OELightGreyPen, 10.0) if args.image: oedepict.OEWriteImage(args.image, image) else: _img = Image.open(io.BytesIO(oedepict.OEWriteImageToBytes("png", image))) _img.show() return os.EX_OK def depict_bfactor_map( image: oedepict.OEImageBase, protein: oechem.OEMolBase, ligand: oechem.OEMolBase, opts: oegrapheme.OE2DActiveSiteDisplayOptions, ) -> None: """Depict B-factor map of active site.""" # perceive interactions active_site = oechem.OEInteractionHintContainer(protein, ligand) if not active_site.IsValid(): oechem.OEThrow.Fatal("Cannot initialize active site!") active_site.SetTitle(ligand.GetTitle()) oechem.OEPerceiveInteractionHints(active_site) # depiction oegrapheme.OEPrepareActiveSiteDepiction(active_site) active_site_disp = oegrapheme.OE2DActiveSiteDisplay(active_site, opts) oegrapheme.OERenderBFactorMap(image, active_site_disp) def get_protein_and_ligand_from_pdb( pdb_filename: str, ) -> tuple[oechem.OEMolBase, oechem.OEMolBase]: """Read protein and and ligand from from pdb/cif file.""" ifs = oechem.oemolistream() if not ifs.open(pdb_filename): oechem.OEThrow.Fatal(f"Unable to open {pdb_filename} for reading") complex_mol = oechem.OEGraphMol() if not oechem.OEReadMolecule(ifs, complex_mol): oechem.OEThrow.Fatal(f"Unable to read complex from {pdb_filename}") if not oechem.OEHasResidues(complex_mol): oechem.OEPerceiveResidues(complex_mol, oechem.OEPreserveResInfo_All) # separate ligand and protein split_opts = oechem.OESplitMolComplexOptions() ligand = oechem.OEGraphMol() protein = oechem.OEGraphMol() water = oechem.OEGraphMol() other = oechem.OEGraphMol() split_opts.SetProteinFilter( oechem.OEOrRoleSet(split_opts.GetProteinFilter(), split_opts.GetWaterFilter()) ) split_opts.SetWaterFilter( oechem.OEMolComplexFilterFactory(oechem.OEMolComplexFilterCategory_Nothing) ) oechem.OESplitMolComplex(ligand, protein, water, other, complex_mol, split_opts) if ligand.NumAtoms() == 0: oechem.OEThrow.Fatal("Cannot separate complex!") return protein, ligand def get_protein_and_ligand_from_design_unit( filename: str, ) -> tuple[oechem.OEMolBase, oechem.OEMolBase]: """Read protein and and ligand from from design unit file.""" du = oechem.OEDesignUnit() if not oechem.OEIsReadableDesignUnit(filename) or not oechem.OEReadDesignUnit( filename, du ): oechem.OEThrow.Fatal("Cannot read design unit.") protein = oechem.OEGraphMol() if not du.GetComponents(protein, oechem.OEDesignUnitComponents_TargetComplex): oechem.OEThrow.Fatal("Could not extract protein from the design unit.") ligand = oechem.OEGraphMol() if not du.GetLigand(ligand): oechem.OEThrow.Fatal("Could not extract ligand from the design unit.") return (protein, ligand) def _check_image_file(args: argparse.Namespace) -> None: # script will terminate if there is some issues if not args.image: # image will be displayed on the screen return ext = Path(args.image).suffix[1:].upper() if not oedepict.OEIsRegisteredImageFile(ext): oechem.OEThrow.Fatal("Unknown image output type!") ofs = oechem.oeofstream() if not ofs.open(args.image): oechem.OEThrow.Fatal("Cannot open output image file!") setattr(main, "__SCRIPT_NAME__", __SCRIPT_NAME__) setattr(main, "__SCRIPT_DESC__", __SCRIPT_DESC__) setattr(main, "__SCRIPT_TOOLKITS__", __SCRIPT_TOOLKITS__) setattr(main, "__SCRIPT_CATEGORIES__", __SCRIPT_CATEGORIES__) if __name__ == "__main__": sys.exit(main())