diff --git a/setup.py b/setup.py index 97d80d313..61ae05cd4 100644 --- a/setup.py +++ b/setup.py @@ -79,7 +79,6 @@ "mrbump = dlstbx.wrapper.mrbump:MrBUMPWrapper", "pandda_xchem = dlstbx.wrapper.pandda_xchem:PanDDAWrapper", "pandda_post = dlstbx.wrapper.pandda_post:PanDDApostWrapper", - "pandda_rhofit = dlstbx.wrapper.pandda_rhofit:PanDDARhofitWrapper", "pipedream_xchem = dlstbx.wrapper.pipedream_xchem:PipedreamWrapper", "phaser_ellg = dlstbx.wrapper.phaser_ellg:PhasereLLGWrapper", "rlv = dlstbx.wrapper.rlv:RLVWrapper", diff --git a/src/dlstbx/util/mvs/helpers.py b/src/dlstbx/util/mvs/helpers.py new file mode 100644 index 000000000..0bddd0a98 --- /dev/null +++ b/src/dlstbx/util/mvs/helpers.py @@ -0,0 +1,55 @@ +from __future__ import annotations + +from pathlib import Path + +import gemmi + + +def find_residue_by_name(structure, name): + for model in structure: + for chain in model: + for res in chain: + if res.name == name: + return chain, res + raise ValueError(f"Residue {name} not found") + + +def residue_centroid(residue): + n = 0 + x = y = z = 0.0 + for at in residue: + p = at.pos + x += p.x + y += p.y + z += p.z + n += 1 + if n == 0: + raise ValueError("Residue has no atoms") + return gemmi.Position(x / n, y / n, z / n) + + +def save_cropped_map(pdb_file, map_file, resname, radius): + st = gemmi.read_structure(pdb_file) + cell = st.cell + m = gemmi.read_ccp4_map(map_file, setup=True) + grid = m.grid + + chain, res = find_residue_by_name(st, resname) + center = residue_centroid(res) # ligand center + + mask = grid.clone() + mask.fill(0.0) + + mask.set_points_around(center, radius, 1.0, use_pbc=True) # spherical mask in Å + + dl = gemmi.Position(radius, radius, radius) # box d/2 + box = gemmi.FractionalBox() + box.extend(cell.fractionalize(center - dl)) + box.extend(cell.fractionalize(center + dl)) + + grid.array[:] *= mask.array + m.set_extent(box) + path = Path(map_file) + map_out = str(path.parents[0] / f"{path.stem}_cropped.ccp4") + m.write_ccp4_map(map_out) + return map_out diff --git a/src/dlstbx/util/mvs/ligandfit.py b/src/dlstbx/util/mvs/viewer_ligandfit.py similarity index 66% rename from src/dlstbx/util/mvs/ligandfit.py rename to src/dlstbx/util/mvs/viewer_ligandfit.py index f5ce216c5..b725e7db7 100644 --- a/src/dlstbx/util/mvs/ligandfit.py +++ b/src/dlstbx/util/mvs/viewer_ligandfit.py @@ -1,10 +1,17 @@ from __future__ import annotations +import gemmi import molviewspec as mvs +from dlstbx.util.mvs.helpers import find_residue_by_name + + +def gen_html_ligandfit(pdb_file, map_file, resname, outdir, acr, smiles, cc): + # make an mvs story from snapshots + st = gemmi.read_structure(pdb_file) + chain, res = find_residue_by_name(st, resname) + residue = mvs.ComponentExpression(label_seq_id=res.seqid.num) -def gen_html_ligandfit(pdb_file, map_file, outdir, acr, smiles, cc): - # make a story from snapshots builder = mvs.create_builder() structure = builder.download(url=pdb_file).parse(format="pdb").model_structure() structure.component(selector="polymer").representation( @@ -30,9 +37,9 @@ def gen_html_ligandfit(pdb_file, map_file, outdir, acr, smiles, cc): snapshot1 = builder.get_snapshot( title="Main View", - description=f"## Ligand_Fit Results: \n ### {acr} with ligand & electron density map \n - SMILES: {smiles} \n - 2FO-FC at 1.5σ, blue \n - Fitting CC = {cc}", - transition_duration_ms=2000, - linger_duration_ms=5000, + description=f"## Ligand_Fit Results: \n ### {acr} with ligand & electron density map \n - SMILES: {smiles} \n - 2FO-FC map at 1.5σ, blue \n - Fitting CC = {cc}", + transition_duration_ms=700, + linger_duration_ms=4000, ) # SNAPSHOT2 @@ -40,11 +47,18 @@ def gen_html_ligandfit(pdb_file, map_file, outdir, acr, smiles, cc): structure = builder.download(url=pdb_file).parse(format="pdb").model_structure() structure.component(selector="polymer").representation( type="surface", size_factor=0.7 - ).opacity(opacity=0.5).color(color="#D8BFD8") - structure.component(selector="polymer").representation().opacity(opacity=0.6).color( - color="grey" + ).opacity(opacity=0.2).color(color="#AABDF1") + structure.component(selector="polymer").representation().opacity( + opacity=0.25 + ).color(custom={"molstar_color_theme_name": "chain_id"}) + structure.component(selector="ligand").representation(type="ball_and_stick").color( + custom={"molstar_color_theme_name": "element-symbol"} ) - structure.component(selector="ligand").focus().representation( + structure.component(selector="ligand").representation(type="surface").opacity( + opacity=0.1 + ).color(custom={"molstar_color_theme_name": "element-symbol"}) + + structure.component(selector=residue).focus().representation( type="ball_and_stick" ).color(custom={"molstar_color_theme_name": "element-symbol"}) @@ -56,25 +70,19 @@ def gen_html_ligandfit(pdb_file, map_file, outdir, acr, smiles, cc): show_faces=False, ).color(color="blue").opacity(opacity=0.25) - # add a label - # info = get_chain_and_residue_numbers(pdb_file, "LIG") - # resid = info[0][1] - residue = mvs.ComponentExpression(label_seq_id=202) - ( - structure.component( - selector=residue, - custom={ - "molstar_show_non_covalent_interactions": True, - "molstar_non_covalent_interactions_radius_ang": 5.0, - }, - ).label(text=f"CC = {cc}") + structure.component( + selector=residue, + custom={ + "molstar_show_non_covalent_interactions": True, + "molstar_non_covalent_interactions_radius_ang": 5, + }, ) snapshot2 = builder.get_snapshot( title="Focus View", description=f"## Ligand_Fit Results: \n ### {acr} with ligand & electron density map \n - SMILES: {smiles} \n - 2FO-FC at 1.5σ, blue \n - Fitting CC = {cc}", - transition_duration_ms=2000, - linger_duration_ms=5000, + transition_duration_ms=700, + linger_duration_ms=4000, ) states = mvs.States( diff --git a/src/dlstbx/util/mvs/viewer_pandda.py b/src/dlstbx/util/mvs/viewer_pandda.py new file mode 100644 index 000000000..0be5d443d --- /dev/null +++ b/src/dlstbx/util/mvs/viewer_pandda.py @@ -0,0 +1,120 @@ +from __future__ import annotations + +from pathlib import Path + +import gemmi +import molviewspec as mvs + +from dlstbx.util.mvs.helpers import find_residue_by_name + + +def gen_html_pandda(pdb_file, event_map, z_map, resname, outdir, dtag, smiles, score): + # make an mvs story from snapshots + st = gemmi.read_structure(pdb_file) + chain, res = find_residue_by_name(st, resname) + residue = mvs.ComponentExpression(label_seq_id=res.seqid.num) + + builder = mvs.create_builder() + structure = builder.download(url=pdb_file).parse(format="pdb").model_structure() + structure.component(selector="polymer").representation( + type="surface", size_factor=0.7 + ).opacity(opacity=0.2).color(color="#AABDF1") + structure.component(selector="polymer").representation().opacity( + opacity=0.25 + ).color(custom={"molstar_color_theme_name": "chain_id"}) + structure.component(selector="ligand").representation( + type="surface", size_factor=0.7 + ).opacity(opacity=0.1).color(custom={"molstar_color_theme_name": "element-symbol"}) + + structure.component(selector=residue).focus().representation( + type="ball_and_stick" + ).color(custom={"molstar_color_theme_name": "element-symbol"}) + + ccp4 = builder.download(url=event_map).parse(format="map") + ccp4.volume().representation( + type="isosurface", + relative_isovalue=3, + show_wireframe=True, + show_faces=False, + ).color(color="blue").opacity(opacity=0.25) + + structure.component( + selector=residue, + custom={ + "molstar_show_non_covalent_interactions": True, + "molstar_non_covalent_interactions_radius_ang": 5, + }, + ) + + snapshot1 = builder.get_snapshot( + title="Event_map", + description=f"## PanDDA2 Results: \n ### {dtag} \n - Ligand score: {score} \n - SMILES: {smiles} \n - Event map at 3σ, blue", + transition_duration_ms=700, + linger_duration_ms=4000, + ) + + # SNAPSHOT2 + builder = mvs.create_builder() + structure = builder.download(url=pdb_file).parse(format="pdb").model_structure() + structure.component(selector="polymer").representation( + type="surface", size_factor=0.7 + ).opacity(opacity=0.2).color(color="#AABDF1") + structure.component(selector="polymer").representation().opacity( + opacity=0.25 + ).color(custom={"molstar_color_theme_name": "chain_id"}) + structure.component(selector="ligand").representation( + type="surface", size_factor=0.7 + ).opacity(opacity=0.1).color(custom={"molstar_color_theme_name": "element-symbol"}) + + structure.component(selector=residue).focus().representation( + type="ball_and_stick" + ).color(custom={"molstar_color_theme_name": "element-symbol"}) + + ccp4 = builder.download(url=z_map).parse(format="map") + ccp4.volume().representation( + type="isosurface", + relative_isovalue=3, + show_wireframe=True, + show_faces=False, + ).color(color="green").opacity(opacity=0.25) + + structure.component( + selector=residue, + custom={ + "molstar_show_non_covalent_interactions": True, + "molstar_non_covalent_interactions_radius_ang": 5, + }, + ) + + snapshot2 = builder.get_snapshot( + title="Z_map", + description=f"## PanDDA2 Results: \n ### {dtag} \n - Ligand score: {score} \n - SMILES: {smiles} \n - Z_map at 3σ, green", + transition_duration_ms=700, + linger_duration_ms=4000, + ) + + states = mvs.States( + snapshots=[snapshot1, snapshot2], # [snapshot1, snapshot2] + metadata=mvs.GlobalMetadata(description="PanDDA2 Results"), + ) + + with open(pdb_file) as f: + pdb_data = f.read() + + # with open(event_map, mode="rb") as f: + # map_data1 = f.read() + + with open(z_map, mode="rb") as f: + map_data2 = f.read() + + html = mvs.molstar_html( + states, + data={pdb_file: pdb_data, z_map: map_data2}, # event_map: map_data1, + ui="stories", + ) + + out_file = Path(f"{outdir}/pandda2_mvs.html") + with open(out_file, "w") as f: + f.write(html) + + return out_file diff --git a/src/dlstbx/util/mvs/viewer_pipedream.py b/src/dlstbx/util/mvs/viewer_pipedream.py new file mode 100644 index 000000000..1473d2743 --- /dev/null +++ b/src/dlstbx/util/mvs/viewer_pipedream.py @@ -0,0 +1,74 @@ +from __future__ import annotations + +from pathlib import Path + +import gemmi +import molviewspec as mvs + +from dlstbx.util.mvs.helpers import find_residue_by_name + + +def gen_html_pipedream(pdb_file, map_file, resname, outdir, dtag, smiles): + # make an mvs story from snapshots + st = gemmi.read_structure(pdb_file) + chain, res = find_residue_by_name(st, resname) + residue = mvs.ComponentExpression(label_seq_id=res.seqid.num) + + builder = mvs.create_builder() + structure = builder.download(url=pdb_file).parse(format="pdb").model_structure() + structure.component(selector="polymer").representation( + type="surface", size_factor=0.7 + ).opacity(opacity=0.2).color(color="#AABDF1") + structure.component(selector="polymer").representation().opacity( + opacity=0.25 + ).color(custom={"molstar_color_theme_name": "chain_id"}) + + structure.component(selector=residue).focus().representation( + type="ball_and_stick" + ).color(custom={"molstar_color_theme_name": "element-symbol"}) + + ccp4 = builder.download(url=map_file).parse(format="map") + ccp4.volume().representation( + type="isosurface", + relative_isovalue=1.5, + show_wireframe=True, + show_faces=False, + ).color(color="blue").opacity(opacity=0.25) + + structure.component( + selector=residue, + custom={ + "molstar_show_non_covalent_interactions": True, + "molstar_non_covalent_interactions_radius_ang": 5, + }, + ) + + snapshot1 = builder.get_snapshot( + title="Main View", + description=f"## Pipedream Results: \n ### {dtag} \n - SMILES: {smiles} \n - 2FO-FC map map at 1.5σ, blue", + transition_duration_ms=700, + linger_duration_ms=4000, + ) + + states = mvs.States( + snapshots=[snapshot1], + metadata=mvs.GlobalMetadata(description="Pipedream Results"), + ) + + with open(pdb_file) as f: + pdb_data = f.read() + + with open(map_file, mode="rb") as f: + map_data = f.read() + + html = mvs.molstar_html( + states, + data={pdb_file: pdb_data, map_file: map_data}, + ui="stories", + ) + + out_file = Path(f"{outdir}/pipedream_mvs.html") + with open(out_file, "w") as f: + f.write(html) + + return out_file diff --git a/src/dlstbx/wrapper/ligand_fit.py b/src/dlstbx/wrapper/ligand_fit.py index a6e3301d5..46b3cee27 100644 --- a/src/dlstbx/wrapper/ligand_fit.py +++ b/src/dlstbx/wrapper/ligand_fit.py @@ -9,7 +9,7 @@ from shutil import ignore_patterns import dlstbx.util.symlink -from dlstbx.util.mvs.ligandfit import gen_html_ligandfit +from dlstbx.util.mvs.viewer_ligandfit import gen_html_ligandfit from dlstbx.wrapper import Wrapper @@ -149,7 +149,13 @@ def run(self): CC = self.pull_CC_from_log(pipeline_directory) try: gen_html_ligandfit( - out_pdb, out_map, pipeline_directory, cc=CC, smiles=smiles, acr=acr + out_pdb, + out_map, + resname="LIG", + outdir=pipeline_directory, + cc=CC, + smiles=smiles, + acr=acr, ) except Exception as e: self.log.debug(f"Exception generating mvs html: {e}") diff --git a/src/dlstbx/wrapper/pandda_rhofit.py b/src/dlstbx/wrapper/pandda_rhofit.py deleted file mode 100644 index 1b23f94fa..000000000 --- a/src/dlstbx/wrapper/pandda_rhofit.py +++ /dev/null @@ -1,301 +0,0 @@ -from __future__ import annotations - -import json -import os -import subprocess -from pathlib import Path - -import gemmi -import numpy as np -import yaml - -from dlstbx.wrapper import Wrapper - - -class PanDDARhofitWrapper(Wrapper): - _logger_name = "dlstbx.wrap.pandda_rhofit" - - def run(self): - assert hasattr(self, "recwrap"), "No recipewrapper object found" - self.log.info( - f"Running recipewrap file {self.recwrap.recipe_step['parameters']['recipewrapper']}" - ) - - params = self.recwrap.recipe_step["job_parameters"] - slurm_task_id = os.environ.get("SLURM_ARRAY_TASK_ID") - # self.log.info((f"SLURM_ARRAY_TASK_ID: {slurm_task_id}")) - - # EVENT_MAP_PATTERN = "{dtag}-event_{event_idx}_1-BDC_{bdc}_map.native.ccp4" - # GROUND_STATE_PATTERN = "{dtag}-ground-state-average-map.native.ccp4" - PANDDA_2_DIR = "/dls_sw/i04-1/software/PanDDA2" - - processing_dir = Path(params.get("processing_directory")) - analysis_dir = processing_dir / "analysis" - model_dir = analysis_dir / "auto_model_building" - auto_panddas_dir = analysis_dir / "auto_pandda2" - - n_datasets = int(params.get("n_datasets")) - self.log.info(f"N_datasets: {n_datasets}") - if n_datasets > 1: - with open(model_dir / ".batch.json", "r") as f: - datasets = json.load(f) - dtag = datasets[int(slurm_task_id) - 1] - else: - dtag = params.get("dtag") - - dataset_dir = auto_panddas_dir / "processed_datasets" / dtag - modelled_dir = dataset_dir / "modelled_structures" - out_dir = modelled_dir / "rhofit" - out_dir.mkdir(parents=True, exist_ok=True) - - self.log.info(f"Processing dtag: {dtag}") - # ------------------------------------------------------- - - event_yaml = dataset_dir / "events.yaml" - - with open(event_yaml, "r") as file: - data = yaml.load(file, Loader=yaml.SafeLoader) - - if not data: - self.log.info( - (f"No events in {event_yaml}, can't continue with PanDDA2 Rhofit") - ) - return False - - # Determine which builds to perform. More than one binder is unlikely and score ranks - # well so build the best scoring event of each dataset. - best_key = max(data, key=lambda k: data[k]["Score"]) - best_entry = data[best_key] - - # event_idx = best_key - # bdc = best_entry["BDC"] - coord = best_entry["Centroid"] - - dataset_dir = auto_panddas_dir / "processed_datasets" / dtag - ligand_dir = dataset_dir / "ligand_files" - build_dmap = dataset_dir / f"{dtag}-z_map.native.ccp4" - restricted_build_dmap = dataset_dir / "build.ccp4" - pdb_file = dataset_dir / f"{dtag}-pandda-input.pdb" - mtz_file = dataset_dir / f"{dtag}-pandda-input.mtz" - restricted_pdb_file = dataset_dir / "build.pdb" - - dmap_cut = 2.0 - # This is usually quite a good contour for building and consistent - # (usually) with the cutoffs PanDDA 2 uses for event finding - - # Rhofit can be confused by hunting non-binding site density. This can be avoided - # by truncating the map to near the binding site - dmap = self.read_pandda_map(build_dmap) - dmap = self.mask_map(dmap, coord) - self.save_xmap(dmap, restricted_build_dmap) - - # Rhofit masks the protein before building. If the original protein - # model clips the event then this results in autobuilding becoming impossible. - # To address tis residues within a 10A neighbourhood of the binding event - # are removed. - self.log.debug("Removing nearby atoms to make room for autobuilding") - self.remove_nearby_atoms( - pdb_file, - coord, - 10.0, - restricted_pdb_file, - ) - - # Really all the cifs should be tried and the best used, or it should try the best - # cif from PanDDA - # This is a temporary fix that will get 90% of situations that can be improved upon - cifs = [x for x in ligand_dir.glob("*.cif")] - if len(cifs) == 0: - self.log.error("No .cif files found!") - - # ------------------------------------------------------- - rhofit_command = f"module load buster; source {PANDDA_2_DIR}/venv/bin/activate; \ - {PANDDA_2_DIR}/scripts/pandda_rhofit.sh -pdb {restricted_pdb_file} -map {build_dmap} -mtz {mtz_file} -cif {cifs[0]} -out {out_dir} -cut {dmap_cut}; " - # cp {modelled_dir}/{dtag}-pandda-model.pdb {modelled_dir}/pandda-internal-fitted.pdb; - - self.log.info(f"Running rhofit command: {rhofit_command}") - - try: - result = subprocess.run( - rhofit_command, - shell=True, - capture_output=True, - text=True, - cwd=auto_panddas_dir, - check=True, - timeout=params.get("timeout-minutes") * 60, - ) - - except subprocess.CalledProcessError as e: - self.log.error(f"Rhofit command: '{rhofit_command}' failed") - self.log.info(e.stdout) - self.log.error(e.stderr) - return False - - with open(out_dir / "rhofit.log", "w") as log_file: - log_file.write(result.stdout) - - # ------------------------------------------------------- - # Merge the protein structure with ligand - protein_st_file = dataset_dir / f"{dtag}-pandda-input.pdb" - ligand_st_file = out_dir / "rhofit" / "best.pdb" - output_file = modelled_dir / f"{dtag}-pandda-model.pdb" - - protein_st = gemmi.read_structure(str(protein_st_file)) - ligand_st = gemmi.read_structure(str(ligand_st_file)) - contact_chain = self.get_contact_chain(protein_st, ligand_st) - protein_st[0][contact_chain].add_residue(ligand_st[0][0][0]) - - protein_st.write_pdb(str(output_file)) - - self.log.info("Auto PanDDA2-Rhofit finished successfully") - return True - - def save_xmap(self, xmap, xmap_file): - """Convenience script for saving ccp4 files.""" - ccp4 = gemmi.Ccp4Map() - ccp4.grid = xmap - ccp4.update_ccp4_header() - ccp4.write_ccp4_map(str(xmap_file)) - - def read_pandda_map(self, xmap_file): - """PanDDA 2 maps are often truncated, and PanDDA 1 maps can have misasigned spacegroups. - This method handles both.""" - dmap_ccp4 = gemmi.read_ccp4_map(str(xmap_file), setup=False) - dmap_ccp4.grid.spacegroup = gemmi.find_spacegroup_by_name("P1") - dmap_ccp4.setup(0.0) - dmap = dmap_ccp4.grid - return dmap - - def expand_event_map(self, bdc, ground_state_file, xmap_file, coord, out_file): - """DEPRECATED. A method for recalculating event maps over the full cell.""" - ground_state_ccp4 = gemmi.read_ccp4_map(str(ground_state_file), setup=False) - ground_state_ccp4.grid.spacegroup = gemmi.find_spacegroup_by_name("P1") - ground_state_ccp4.setup(0.0) - ground_state = ground_state_ccp4.grid - - xmap_ccp4 = gemmi.read_ccp4_map(str(xmap_file), setup=False) - xmap_ccp4.grid.spacegroup = gemmi.find_spacegroup_by_name("P1") - xmap_ccp4.setup(0.0) - xmap = xmap_ccp4.grid - - mask = gemmi.FloatGrid(xmap.nu, xmap.nv, xmap.nw) - mask.set_unit_cell(xmap.unit_cell) - mask.set_points_around( - gemmi.Position(coord[0], coord[1], coord[2]), radius=10.0, value=1.0 - ) - - event_map = gemmi.FloatGrid(xmap.nu, xmap.nv, xmap.nw) - event_map.set_unit_cell(xmap.unit_cell) - event_map_array = np.array(event_map, copy=False) - event_map_array[:, :, :] = np.array(xmap)[:, :, :] - ( - bdc * np.array(ground_state)[:, :, :] - ) - event_map_array[:, :, :] = event_map_array[:, :, :] * np.array(mask)[:, :, :] - - event_map_non_zero = event_map_array[event_map_array != 0.0] - cut = np.std(event_map_non_zero) - - return cut - - def mask_map(self, dmap, coord, radius=10.0): - """Simple routine to mask density to region around a specified point.""" - mask = gemmi.FloatGrid(dmap.nu, dmap.nv, dmap.nw) - mask.set_unit_cell(dmap.unit_cell) - mask.set_points_around( - gemmi.Position(coord[0], coord[1], coord[2]), radius=radius, value=1.0 - ) - - dmap_array = np.array(dmap, copy=False) - dmap_array[:, :, :] = dmap_array[:, :, :] * np.array(mask)[:, :, :] - - return dmap - - def remove_nearby_atoms(self, pdb_file, coord, radius, output_file): - """An inelegant method for removing residues near the event centroid and creating - a new, truncated pdb file. GEMMI doesn't have a super nice way to remove - residues according to a specific criteria.""" - st = gemmi.read_structure(str(pdb_file)) - new_st = st.clone() # Clone to keep metadata - - coord_array = np.array([coord[0], coord[1], coord[2]]) - - # Delete all residues for a clean chain. Yes this is an arcane way to do it. - chains_to_delete = [] - for model in st: - for chain in model: - chains_to_delete.append((model.num, chain.name)) - - for model in new_st: - for chain in model: - for res in chain: - del chain[-1] - - # Add non-rejected residues to a new structure - for j, model in enumerate(st): - for k, chain in enumerate(model): - for res in chain: - add_res = True - for atom in res: - pos = atom.pos - distance = np.linalg.norm( - coord_array - np.array([pos.x, pos.y, pos.z]) - ) - if distance < radius: - add_res = False - - if add_res: - new_st[j][k].add_residue(res) - new_st.write_pdb(str(output_file)) - - def get_contact_chain(self, protein_st, ligand_st): - """A simple estimation of the contact chain based on which chain has the most atoms - nearby.""" - ligand_pos_list = [] - for model in protein_st: - for chain in model: - for res in chain: - for atom in res: - pos = atom.pos - ligand_pos_list.append([pos.x, pos.y, pos.z]) - centroid = np.linalg.norm(np.array(ligand_pos_list), axis=0) - - PROTEIN_RESIDUES = [ - "ALA", - "ARG", - "ASN", - "ASP", - "CYS", - "GLN", - "GLU", - "HIS", - "ILE", - "LEU", - "LYS", - "MET", - "PHE", - "PRO", - "SER", - "THR", - "TRP", - "TYR", - "VAL", - "GLY", - ] - - chain_counts = {} - for model in protein_st: - for chain in model: - chain_counts[chain.name] = 0 - for res in chain: - if res.name not in PROTEIN_RESIDUES: - continue - for atom in res: - pos = atom.pos - distance = np.linalg.norm( - np.array([pos.x, pos.y, pos.z]) - centroid - ) - if distance < 5.0: - chain_counts[chain.name] += 1 - - return min(chain_counts, key=lambda _x: chain_counts[_x]) diff --git a/src/dlstbx/wrapper/pandda_xchem.py b/src/dlstbx/wrapper/pandda_xchem.py index 8865ba4b5..5b20447c8 100644 --- a/src/dlstbx/wrapper/pandda_xchem.py +++ b/src/dlstbx/wrapper/pandda_xchem.py @@ -3,7 +3,6 @@ import json import os import shutil -import sqlite3 import subprocess from pathlib import Path @@ -11,6 +10,12 @@ import numpy as np import yaml +import dlstbx.util.symlink +from dlstbx.util.mvs.helpers import ( + find_residue_by_name, + save_cropped_map, +) +from dlstbx.util.mvs.viewer_pandda import gen_html_pandda from dlstbx.wrapper import Wrapper @@ -32,8 +37,8 @@ def run(self): analysis_dir = Path(processed_dir / "analysis") pandda_dir = analysis_dir / "pandda2" model_dir = pandda_dir / "model_building" - auto_panddas_dir = Path(pandda_dir / "panddas") - Path(auto_panddas_dir).mkdir(exist_ok=True) + panddas_dir = Path(pandda_dir / "panddas") + Path(panddas_dir).mkdir(exist_ok=True) n_datasets = int(params.get("n_datasets")) if n_datasets > 1: # array job case @@ -43,10 +48,20 @@ def run(self): else: dtag = params.get("dtag") - self.log.info(f"Processing dtag: {dtag}") dataset_dir = model_dir / dtag compound_dir = dataset_dir / "compound" + if pipeline_final_params := params.get("pipeline-final", []): + final_directory = Path(pipeline_final_params["path"]) + final_directory.mkdir(parents=True, exist_ok=True) + if params.get("create_symlink"): + dlstbx.util.symlink.create_parent_symlink( + final_directory, params.get("create_symlink") + ) + + self.log.info(f"Processing dtag: {dtag}") + + smiles = params.get("smiles") smiles_files = list(compound_dir.glob("*.smiles")) if len(smiles_files) == 0: @@ -56,7 +71,7 @@ def run(self): return False elif len(smiles_files) > 1: self.log.error( - f"Multiple .smiles files found in in {compound_dir}:, {smiles_files}, warning for dtag {dtag}" + f"Multiple .smiles files found in in {compound_dir}: {smiles_files}, warning for dtag {dtag}" ) return False @@ -64,18 +79,22 @@ def run(self): CompoundCode = smiles_file.stem # ------------------------------------------------------- + # Ligand restraint generation + + restraints_log = dataset_dir / "restraints.log" + attachments = [restraints_log] # synchweb attachments + restraints_command = f"grade2 --in {smiles_file} --itype smi --out {CompoundCode} -f > {restraints_log}" # acedrg_command = f"module load ccp4; acedrg -i {smiles_file} -o {CompoundCode}" - restraints_command = f"grade2 --in {smiles_file} --itype smi --out {CompoundCode} -f" try: - result = subprocess.run( + subprocess.run( restraints_command, shell=True, capture_output=True, text=True, cwd=compound_dir, check=True, - timeout=params.get("timeout-minutes") * 60, + timeout=30 * 60, ) except subprocess.CalledProcessError as e: @@ -85,23 +104,33 @@ def run(self): self.log.info(e.stdout) self.log.error(e.stderr) + self.send_attachments_to_ispyb(attachments, final_directory) return False restraints = compound_dir / f"{CompoundCode}.restraints.cif" restraints.rename(compound_dir / f"{CompoundCode}.cif") - pdb = compound_dir / f"{CompoundCode}.xyz.pdb" - pdb.rename(compound_dir / f"{CompoundCode}.pdb") + ligand_pdb = compound_dir / f"{CompoundCode}.xyz.pdb" + ligand_pdb.rename(compound_dir / f"{CompoundCode}.pdb") - with open(dataset_dir / "restraints.log", "w") as log_file: - log_file.write(result.stdout) + ligand_cif = compound_dir / f"{CompoundCode}.cif" - self.log.info(f"Restraints generated succesfully for dtag {dtag}, launching PanDDA2") + self.log.info( + f"Restraints generated succesfully for dtag {dtag}, launching PanDDA2" + ) + + # ------------------------------------------------------- + # PanDDA2 - pandda2_command = f"source /dls_sw/i04-1/software/PanDDA2/venv/bin/activate; \ - python -u /dls_sw/i04-1/software/PanDDA2/scripts/process_dataset.py --data_dirs={model_dir} --out_dir={auto_panddas_dir} --dtag={dtag} --use_ligand_data=False --local_cpus=1" + dataset_pdir = panddas_dir / "processed_datasets" / dtag + dataset_pdir.mkdir(exist_ok=True) + pandda2_log = dataset_pdir / "pandda2.log" + + attachments.extend([pandda2_log, ligand_cif]) + pandda2_command = f"source {PANDDA_2_DIR}/venv/bin/activate; \ + python -u /dls_sw/i04-1/software/PanDDA2/scripts/process_dataset.py --data_dirs={model_dir} --out_dir={panddas_dir} --dtag={dtag} --use_ligand_data=False --local_cpus=1 > {pandda2_log}" try: - result = subprocess.run( + subprocess.run( pandda2_command, shell=True, capture_output=True, @@ -115,9 +144,12 @@ def run(self): self.log.error(f"PanDDA2 command: '{pandda2_command}' failed") self.log.info(e.stdout) self.log.error(e.stderr) + self.send_attachments_to_ispyb(attachments, final_directory) return False - dataset_pdir = auto_panddas_dir / "processed_datasets" / dtag + # ------------------------------------------------------- + # PanDDA Rhofit ligand fitting + ligand_dir = dataset_pdir / "ligand_files" # pandda2 not moving files into ligand_dir, fix @@ -129,52 +161,47 @@ def run(self): target.unlink() target.symlink_to(file) - pandda_log = dataset_pdir / "pandda2.log" - with open(pandda_log, "w") as log_file: - log_file.write(result.stdout) - modelled_dir = dataset_pdir / "modelled_structures" out_dir = modelled_dir / "rhofit" out_dir.mkdir(parents=True, exist_ok=True) - event_yaml = dataset_pdir / "events.yaml" - - with open(event_yaml, "r") as file: - data = yaml.load(file, Loader=yaml.SafeLoader) + events_yaml = dataset_pdir / "events.yaml" - if not data: + if not events_yaml.exists(): self.log.info( - (f"No events in {event_yaml}, can't continue with PanDDA2 Rhofit") + (f"No events in {events_yaml}, can't continue with PanDDA2 Rhofit") ) - return True # False + self.send_attachments_to_ispyb(attachments, final_directory) + return False + + with open(events_yaml, "r") as file: + data = yaml.load(file, Loader=yaml.SafeLoader) # Determine which builds to perform. More than one binder is unlikely and score ranks # well so build the best scoring event of each dataset. best_key = max(data, key=lambda k: data[k]["Score"]) best_entry = data[best_key] - # event_idx = best_key + event_idx = best_key # bdc = best_entry["BDC"] coord = best_entry["Centroid"] - build_dmap = dataset_pdir / f"{dtag}-z_map.native.ccp4" restricted_build_dmap = dataset_pdir / "build.ccp4" + z_map = dataset_pdir / f"{dtag}-z_map.native.ccp4" + event_map = list(dataset_pdir.glob(f"{dtag}-event_{event_idx}*"))[0] + # event_map = dataset_dir / f'{dtag}-event_{event_idx}_1-BDC_{bdc}_map.native.ccp4' #round BDC to 2dp? pdb_file = dataset_pdir / f"{dtag}-pandda-input.pdb" mtz_file = dataset_pdir / f"{dtag}-pandda-input.mtz" restricted_pdb_file = dataset_pdir / "build.pdb" - dmap_cut = 2.0 - # This is usually quite a good contour for building and consistent - # (usually) with the cutoffs PanDDA 2 uses for event finding - # Rhofit can be confused by hunting non-binding site density. This can be avoided # by truncating the map to near the binding site - dmap = self.read_pandda_map(build_dmap) + dmap = self.read_pandda_map(event_map) dmap = self.mask_map(dmap, coord) self.save_xmap(dmap, restricted_build_dmap) # Rhofit masks the protein before building. If the original protein # model clips the event then this results in autobuilding becoming impossible. - # To address tis residues within a 10A neighbourhood of the binding event + # To address this residues within a 10A neighbourhood of the binding event # are removed. self.log.debug("Removing nearby atoms to make room for autobuilding") self.remove_nearby_atoms( @@ -184,57 +211,116 @@ def run(self): restricted_pdb_file, ) - # Really all the cifs should be tried and the best used, or it should try the best - # cif from PanDDA - # This is a temporary fix that will get 90% of situations that can be improved upon cifs = list(ligand_dir.glob("*.cif")) - if len(cifs) == 0: - self.log.error( - f"No .cif files found for dtag {dtag}, cannot launch PanDDA2 Rhofit!" - ) - return True - # ------------------------------------------------------- + rhofit_log = dataset_pdir / "rhofit.log" + attachments.extend([event_map, z_map, rhofit_log]) rhofit_command = f"module load buster; source {PANDDA_2_DIR}/venv/bin/activate; \ - {PANDDA_2_DIR}/scripts/pandda_rhofit.sh -pdb {restricted_pdb_file} -map {build_dmap} -mtz {mtz_file} -cif {cifs[0]} -out {out_dir} -cut {dmap_cut}; " + {PANDDA_2_DIR}/scripts/pandda_rhofit.sh -pdb {restricted_pdb_file} -map {restricted_build_dmap} -mtz {mtz_file} -cif {cifs[0]} -out {out_dir} -cut 2.0 > {rhofit_log};" # dmap_cut=2.0 self.log.info(f"Running PanDDA Rhofit command: {rhofit_command}") try: - result = subprocess.run( + subprocess.run( rhofit_command, shell=True, capture_output=True, text=True, - cwd=auto_panddas_dir, + cwd=panddas_dir, check=True, - timeout=params.get("timeout-minutes") * 60, + timeout=60 * 60, ) except subprocess.CalledProcessError as e: self.log.error(f"Rhofit command: '{rhofit_command}' failed") self.log.info(e.stdout) self.log.error(e.stderr) + self.send_attachments_to_ispyb(attachments, final_directory) return False - with open(out_dir / "rhofit.log", "w") as log_file: - log_file.write(result.stdout) - - # ------------------------------------------------------- - # Merge the protein structure with ligand + # Merge the protein structure with ligand -> pandda model protein_st_file = dataset_pdir / f"{dtag}-pandda-input.pdb" ligand_st_file = out_dir / "rhofit" / "best.pdb" - output_file = modelled_dir / f"{dtag}-pandda-model.pdb" + pandda_model = modelled_dir / f"{dtag}-pandda-model.pdb" protein_st = gemmi.read_structure(str(protein_st_file)) ligand_st = gemmi.read_structure(str(ligand_st_file)) contact_chain = self.get_contact_chain(protein_st, ligand_st) protein_st[0][contact_chain].add_residue(ligand_st[0][0][0]) - if output_file.exists(): - shutil.copy(output_file, modelled_dir / "pandda-internal-fitted.pdb") + if pandda_model.exists(): # backup previous model + shutil.copy2(pandda_model, modelled_dir / "pandda-internal-fitted.pdb") - protein_st.write_pdb(str(output_file)) + protein_st.write_pdb(str(pandda_model)) + + # ------------------------------------------------------- + # Ligand scoring + + ligand_score = dataset_pdir / "ligand_score.txt" + attachments.extend([pandda_model, ligand_score]) + + st = gemmi.read_structure(str(pandda_model)) + chain, res = find_residue_by_name(st, "LIG") + ligand_id = chain.name + f"/{res.seqid.num}" + + score_command = f"source {PANDDA_2_DIR}/venv/bin/activate; \ + python {PANDDA_2_DIR}/scripts/ligand_score.py --mtz_path={mtz_file} --zmap_path={z_map} --ligand_id={ligand_id} --structure_path={pandda_model} --out_path={ligand_score}" + + self.log.info(f"Running Ligand Score command: {score_command}") + + try: + subprocess.run( + score_command, + shell=True, + capture_output=True, + text=True, + cwd=panddas_dir, + check=True, + timeout=10 * 60, + ) + + except subprocess.CalledProcessError as e: + self.log.error(f"Ligand score command: '{score_command}' failed") + self.log.info(e.stdout) + self.log.error(e.stderr) + + with open(ligand_score, "r") as file: + score = float(file.read().strip()) + + self.log.info(f"Ligand score for {dtag} = {score}") + + try: + cropped_event_map = save_cropped_map( + str(pandda_model), str(z_map), "LIG", radius=6 + ) + cropped_z_map = save_cropped_map( + str(pandda_model), str(z_map), "LIG", radius=6 + ) + mvs_html = gen_html_pandda( + str(pandda_model), + cropped_event_map, + cropped_z_map, + resname="LIG", + outdir=dataset_pdir, + dtag=dtag, + smiles=smiles, + score=score, + ) + attachments.extend([mvs_html]) + except Exception as e: + self.log.debug(f"Exception generating mvs html: {e}") + + # data = [ + # ["SMILES code", "Fitted_ligand?", "Hit?"], + # [f"{smiles}", "True", "True"], + # ] + # json_results = dataset_pdir / "pandda2_results.json" + # with open(json_results, "w") as f: + # json.dump(data, f) + # attachments.extend([json_results]) + + self.log.info(f"Attachments list: {attachments}") + self.send_attachments_to_ispyb(attachments, final_directory) self.log.info(f"Auto PanDDA2 pipeline finished successfully for dtag {dtag}") return True @@ -299,7 +385,7 @@ def mask_map(self, dmap, coord, radius=10.0): return dmap - def remove_nearby_atoms(self, pdb_file, coord, radius, output_file): + def remove_nearby_atoms(self, pdb_file, coord, radius, pandda_model): """An inelegant method for removing residues near the event centroid and creating a new, truncated pdb file. GEMMI doesn't have a super nice way to remove residues according to a specific criteria.""" @@ -334,7 +420,7 @@ def remove_nearby_atoms(self, pdb_file, coord, radius, output_file): if add_res: new_st[j][k].add_residue(res) - new_st.write_pdb(str(output_file)) + new_st.write_pdb(str(pandda_model)) def get_contact_chain(self, protein_st, ligand_st): """A simple estimation of the contact chain based on which chain has the most atoms @@ -388,14 +474,38 @@ def get_contact_chain(self, protein_st, ligand_st): return min(chain_counts, key=lambda _x: chain_counts[_x]) - def update_data_source(self, db_dict, dtag, database_path): - sql = ( - "UPDATE mainTable SET " - + ", ".join([f"{k} = :{k}" for k in db_dict]) - + f" WHERE CrystalName = '{dtag}'" - ) - conn = sqlite3.connect(database_path) - # conn.execute("PRAGMA journal_mode=WAL;") - cursor = conn.cursor() - cursor.execute(sql, db_dict) - conn.commit() + def send_attachments_to_ispyb(self, attachments, final_directory): + for f in attachments: + if f.exists(): + if f.suffix == ".html": + file_type = "Result" + importance_rank = 1 + elif f.suffix == ".ccp4": + file_type = "Result" + importance_rank = 1 + elif f.suffix == ".cif": + file_type = "Result" + importance_rank = 1 + elif f.suffix == ".pdb": + file_type = "Result" + importance_rank = 1 + elif f.suffix == ".log": + file_type = "Log" + importance_rank = 2 + else: + continue + try: + shutil.copy(f, final_directory) + result_dict = { + "file_path": str(final_directory), + "file_name": f.name, + "file_type": file_type, + "importance_rank": importance_rank, + } + self.record_result_individual_file(result_dict) + self.log.info(f"Uploaded {f.name} as an attachment") + + except Exception: + self.log.warning( + f"Could not attach {f.name} to ISPyB", exc_info=True + ) diff --git a/src/dlstbx/wrapper/pipedream_xchem.py b/src/dlstbx/wrapper/pipedream_xchem.py index 13f27a3db..1eeda8dee 100644 --- a/src/dlstbx/wrapper/pipedream_xchem.py +++ b/src/dlstbx/wrapper/pipedream_xchem.py @@ -2,12 +2,16 @@ import json import os +import shutil import sqlite3 import subprocess from pathlib import Path import portalocker +import dlstbx.util.symlink +from dlstbx.util.mvs.helpers import save_cropped_map +from dlstbx.util.mvs.viewer_pipedream import gen_html_pipedream from dlstbx.wrapper import Wrapper @@ -37,6 +41,14 @@ def run(self): dimple_mtz = dataset_dir / "dimple.mtz" upstream_mtz = dataset_dir / f"{dtag}.free.mtz" + if pipeline_final_params := params.get("pipeline-final", []): + final_directory = Path(pipeline_final_params["path"]) + final_directory.mkdir(parents=True, exist_ok=True) + if params.get("create_symlink"): + dlstbx.util.symlink.create_parent_symlink( + final_directory, params.get("create_symlink") + ) + self.log.info(f"Processing dtag: {dtag}") dataset_dir = model_dir / dtag @@ -59,39 +71,39 @@ def run(self): CompoundCode = smiles_file.stem # ------------------------------------------------------- - restraints_command = f"module load buster; module load graphviz; \ - export CSDHOME=/dls_sw/apps/CSDS/2024.1.0/; export BDG_TOOL_MOGUL=/dls_sw/apps/CSDS/2024.1.0/ccdc-software/mogul/bin/mogul; \ - grade2 --in {smiles_file} --itype smi --out {CompoundCode} -f" + # Ligand restraint generation + + restraints_log = dataset_dir / "restraints.log" + attachments = [restraints_log] # synchweb attachments + + restraints_command = f"grade2 --in {smiles_file} --itype smi --out {CompoundCode} -f > {restraints_log}" try: - result = subprocess.run( + subprocess.run( restraints_command, shell=True, capture_output=True, text=True, cwd=compound_dir, check=True, - timeout=params.get("timeout-minutes") * 60, + timeout=30 * 60, ) except subprocess.CalledProcessError as e: self.log.error( f"Ligand restraint generation command: '{restraints_command}' failed for dataset {dtag}" ) - self.log.info(e.stdout) self.log.error(e.stderr) + self.send_attachments_to_ispyb(attachments, final_directory) return False restraints = compound_dir / f"{CompoundCode}.restraints.cif" restraints.rename(compound_dir / f"{CompoundCode}.cif") - pdb = compound_dir / f"{CompoundCode}.xyz.pdb" - pdb.rename(compound_dir / f"{CompoundCode}.pdb") - - with open(dataset_dir / "restraints.log", "w") as log_file: - log_file.write(result.stdout) + ligand_pdb = compound_dir / f"{CompoundCode}.xyz.pdb" + ligand_pdb.rename(compound_dir / f"{CompoundCode}.pdb") - ligand_cif = str(compound_dir / f"{CompoundCode}.cif") + ligand_cif = compound_dir / f"{CompoundCode}.cif" self.log.info(f"Restraints generated succesfully for dtag {dtag}") self.log.info(f"Removing crystallisation components from pdb file for {dtag}") @@ -99,10 +111,12 @@ def run(self): self.log.info(f"Launching pipedream for {dtag}") # ------------------------------------------------------- + # Pipedream - pipedream_command = f"module load buster; module load graphviz; \ - export BDG_TOOL_MOGUL=/dls_sw/apps/CSDS/2024.1.0/ccdc-software/mogul/bin/mogul; \ - /dls_sw/apps/GPhL/BUSTER/20250717/scripts/pipedream \ + pipedream_log = out_dir / "summary.out" + attachments.extend([pipedream_log, ligand_cif]) + + pipedream_command = f"/dls_sw/apps/GPhL/BUSTER/20250717/scripts/pipedream \ -nolmr \ -hklin {upstream_mtz} \ -xyzin {dimple_pdb} \ @@ -116,10 +130,10 @@ def run(self): -rhocommands \ -xclusters \ -nochirals \ - -rhofit {ligand_cif}" + -rhofit {str(ligand_cif)}" try: - result = subprocess.run( + subprocess.run( pipedream_command, shell=True, capture_output=True, @@ -133,11 +147,26 @@ def run(self): self.log.error(f"Pipedream command: '{pipedream_command}' failed") self.log.info(e.stdout) self.log.error(e.stderr) + + with open(out_dir / "stderr.out", "w") as stderr: + stderr.write(e.stderr) + + attachments.extend([out_dir / "stderr.out"]) + self.send_attachments_to_ispyb(attachments, final_directory) return False self.log.info(f"Pipedream finished successfully for dtag {dtag}") - pipedream_summary = f"{out_dir}/pipedream_summary.json" + # ------------------------------------------------------- + + report_dir = out_dir / f"report-{CompoundCode}" + buster_report = report_dir / "report.pdf" + + postrefine_dir = out_dir / f"postrefine-{CompoundCode}" + refine_mtz = postrefine_dir / "refine.mtz" + refine_pdb = postrefine_dir / "refine.pdb" + + pipedream_summary = out_dir / "pipedream_summary.json" self.save_dataset_metadata( str(pipedream_dir), str(compound_dir), @@ -148,6 +177,8 @@ def run(self): dtag, ) + attachments.extend([buster_report, refine_mtz, refine_pdb, pipedream_summary]) + try: with open(pipedream_summary, "r") as f: data = json.load(f) @@ -163,12 +194,10 @@ def run(self): ) except Exception as e: self.log.info(f"Can't continue with pipedream postprocessing: {e}") + self.send_attachments_to_ispyb(attachments, final_directory) return True # Post-processing: Generate maps and run edstats - postrefine_dir = out_dir / f"postrefine-{CompoundCode}" - refine_mtz = postrefine_dir / "refine.mtz" - refine_pdb = postrefine_dir / "refine.pdb" map_2fofc = postrefine_dir / "refine_2fofc.map" map_fofc = postrefine_dir / "refine_fofc.map" @@ -177,12 +206,30 @@ def run(self): os.system(f"gemmi sf2map --sample 5 {str(refine_mtz)} {map_fofc} 2>&1") except Exception as e: self.log.debug(f"Cannot continue with pipedream postprocessing: {e}") + self.send_attachments_to_ispyb(attachments, final_directory) return True + try: + cropped_map = save_cropped_map( + str(refine_pdb), str(map_2fofc), "LIG", radius=6 + ) + mvs_html = gen_html_pipedream( + str(refine_pdb), + cropped_map, + resname="LIG", + outdir=out_dir, + dtag=dtag, + smiles=smiles, + ) + attachments.extend([mvs_html]) + except Exception as e: + self.log.debug(f"Exception generating mvs html: {e}") + if reslo is None or reshi is None: self.log.debug( "Can't continue with pipedream postprocessing: resolution range None" ) + self.send_attachments_to_ispyb(attachments, final_directory) return True # Run edstats if both maps exist and resolution range is found @@ -190,13 +237,15 @@ def run(self): self.log.debug( "Can't continue with pipedream postprocessing: maps not found" ) + self.send_attachments_to_ispyb(attachments, final_directory) return True - edstats_command = f"module load ccp4; edstats XYZIN {refine_pdb} MAPIN1 {map_2fofc} MAPIN2 {map_fofc} OUT {str(postrefine_dir / 'edstats.out')}" + edstats_out = postrefine_dir / "edstats.out" + edstats_command = f"module load ccp4; edstats XYZIN {refine_pdb} MAPIN1 {map_2fofc} MAPIN2 {map_fofc} OUT {str(edstats_out)}" stdin_text = f"RESLO={reslo}\nRESHI={reshi}\nEND\n" try: - result = subprocess.run( + subprocess.run( edstats_command, input=stdin_text, text=True, @@ -209,9 +258,13 @@ def run(self): self.log.error(f"Edstats command: '{edstats_command}' failed") self.log.info(e.stdout) self.log.error(e.stderr) + self.send_attachments_to_ispyb(attachments, final_directory) return True self.log.info(f"Pipedream postprocessing finished successfully for dtag {dtag}") + + attachments.extend([edstats_out]) + self.send_attachments_to_ispyb(attachments, final_directory) return True def process_pdb_file(self, dimple_pdb: Path): @@ -299,38 +352,47 @@ def save_dataset_metadata( with open(json_file, "w", encoding="utf-8") as f: json.dump(data, f, indent=2) - # def send_attachments_to_ispyb(self, pipeline_directory, final_directory): - # for f in pipeline_directory.iterdir(): - # if f.stem.endswith("final"): - # file_type = "Result" - # importance_rank = 1 - # elif f.suffix == ".html": - # file_type = "Result" - # importance_rank = 1 - # elif f.suffix == ".png": - # file_type = "Result" - # importance_rank = 1 - # elif f.suffix == ".json": - # file_type = "Result" - # importance_rank = 1 - # elif f.suffix == ".log": - # file_type = "Log" - # importance_rank = 2 - # else: - # continue - # try: - # shutil.copy(pipeline_directory / f.name, final_directory) - # result_dict = { - # "file_path": str(final_directory), - # "file_name": f.name, - # "file_type": file_type, - # "importance_rank": importance_rank, - # } - # self.record_result_individual_file(result_dict) - # self.log.info(f"Uploaded {f.name} as an attachment") - - # except Exception: - # self.log.warning(f"Could not attach {f.name} to ISPyB", exc_info=True) + def send_attachments_to_ispyb(self, attachments, final_directory): + for f in attachments: + if f.exists(): + if f.suffix == ".html": + file_type = "Result" + importance_rank = 1 + elif f.suffix == ".mtz": + file_type = "Result" + importance_rank = 1 + elif f.suffix == ".cif": + file_type = "Result" + importance_rank = 1 + elif f.suffix == ".pdb": + file_type = "Result" + importance_rank = 1 + elif f.suffix == ".pdf": + file_type = "Result" + importance_rank = 1 + elif f.suffix == ".out": + file_type = "Log" + importance_rank = 2 + elif f.suffix == ".log": + file_type = "Log" + importance_rank = 2 + else: + continue + try: + shutil.copy(f, final_directory) + result_dict = { + "file_path": str(final_directory), + "file_name": f.name, + "file_type": file_type, + "importance_rank": importance_rank, + } + self.record_result_individual_file(result_dict) + self.log.info(f"Uploaded {f.name} as an attachment") + + except Exception: + self.log.warning( + f"Could not attach {f.name} to ISPyB", exc_info=True + ) def update_data_source(self, db_dict, dtag, database_path): sql = (