Coot Scripting Interface  7000
Classes | Functions
cc-interface.hh File Reference

Coot Scripting Interface - General (C++ functions) More...

Go to the source code of this file.

Classes

class  coot::alias_path_t
 alias path More...
 
class  coot::pisa_interface_bond_info_t
 pisa internal function More...
 
class  coot::str_mtime
 str mtime for for attributes More...
 
class  coot::file_attribs_info_t
 trivial helper function for file attributes More...
 

Functions

int rigid_body_fit_with_residue_ranges (int imol, const std::vector< coot::residue_range_t > &ranges)
 return 0 on fail to refine (no sensible place to put atoms) and 1 on fitting happened.
 
int morph_fit_chain (int imol, std::string chain_id, float transformation_averaging_radius)
 Morph the given chain.
 
int morph_fit_residues (int imol, const std::vector< coot::residue_spec_t > &residue_specs, float transformation_averaging_radius)
 morph the given residues.
 
int morph_fit_by_secondary_structure_elements (int imol, const std::string &chain_id)
 morph transformation are based primarily on rigid body refinement of the secondary structure elements.
 
SCM find_blobs_scm (int imol_model, int imol_map, float cut_off_density_level)
 find blobs
 
void b_factor_distribution_graph (int imol)
 B-factor distribution histogram.
 
void get_coords_for_accession_code (const std::string &code)
 if possible, read in the new coords getting coords via web. More...
 
void set_radial_map_colouring_centre (int imol, float x, float y, float z)
 radial map colouring centre
 
void set_radial_map_colouring_min_radius (int imol, float r)
 radial map colouring min
 
void set_radial_map_colouring_max_radius (int imol, float r)
 radial map colouring max
 
void set_radial_map_colouring_invert (int imol, int invert_state)
 radial map colouring inverted colour map
 
void set_radial_map_colouring_saturation (int imol, float saturation)
 radial map colouring saturation More...
 
void prodrg_import_function (std::string file_name, std::string comp_id)
 import given mdl file into prodrg or other 3d generation program More...
 
void sbase_import_function (std::string comp_id)
 import molecule from CCP4 SRS (or SBase, as it used to be called). More...
 
std::pair< int, std::string > align_to_closest_chain (std::string target_seq, float match_fraction)
 align sequence to closest chain (compare across all chains in all molecules). More...
 
void simple_text_dialog (const std::string &dialog_title, const std::string &text, int geom_x, int geom_y)
 make a simple text dialog.
 
void graphics_to_phenix_geo_representation (int imol, int mode, const coot::phenix_geo_bonds &g)
 phenix GEO bonds representation
 
void graphics_to_phenix_geo_representation (int imol, int mode, const std::string &geo_file_name)
 phenix GEO bonds representation, read from file
 
void set_python_draw_function (const std::string &command_string)
 client/server functions
 
int encode_ints (int i1, int i2)
 encoding of ints
 
void store_keyed_user_name (std::string key, std::string user_name, std::string passwd)
 store username and password for the database.
 
More Symmetry Functions
*SCM get_symmetry (int imol)
 return the symmetry of the imolth molecule More...
 
PyObject * get_symmetry_py (int imol)
 
int clashes_with_symmetry (int imol, const char *chain_id, int res_no, const char *ins_code, float clash_dist)
 return 1 if this residue clashes with the symmetry-related atoms of the same molecule. More...
 
void add_molecular_symmetry (int imol, double r_00, double r_01, double r_02, double r_10, double r_11, double r_12, double r_20, double r_21, double r_22, double about_origin_x, double about_origin_y, double about_origin_z)
 
int add_molecular_symmetry_from_mtrix_from_file (int imol, const std::string &file_name)
 
int add_molecular_symmetry_from_mtrix_from_self_file (int imol)
 
Extra Map Functions
std::vector< int > auto_read_make_and_draw_maps (const char *filename)
 read MTZ file filename and from it try to make maps More...
 
std::vector< int > auto_read_make_and_draw_maps_from_mtz (const char *filename)
 set the flag to do a difference map (too) on auto-read MTZ
 
std::vector< int > auto_read_make_and_draw_maps_from_cns (const char *filename)
 
void add_map_colour_mol_menu_item (int imol, const std::string &name, GtkWidget *sub_menu, GtkSignalFunc callback)
 
void add_map_scroll_wheel_mol_menu_item (int imol, const std::string &name, GtkWidget *menu, GtkSignalFunc callback)
 
int sharpen_blur_map (int imol_map, float b_factor)
 make a sharpened or blurred map More...
 
int sharpen_blur_map_with_resampling (int imol_map, float b_factor, float resample_factor)
 make a sharpened or blurred map with resampling More...
 
void multi_sharpen_blur_map_scm (int imol_map, SCM b_factors_list)
 make many sharpened or blurred maps More...
 
void multi_sharpen_blur_map_py (int imol_map, PyObject *b_factors_list)
 make many sharpened or blurred maps More...
 
PyObject * amplitude_vs_resolution_py (int mol_map)
 
SCM amplitude_vs_resolution_scm (int mol_map)
 
int flip_hand (int imol_map)
 Flip the hand of the map. More...
 
int analyse_map_point_density_change (const std::vector< int > &map_number_list)
 test function for analysis of multiple map
 
int analyse_map_point_density_change_py (PyObject *map_number_list)
 
void go_to_map_molecule_centre (int imol_map)
 Go to the centre of the molecule - for Cryo-EM Molecules. More...
 
float b_factor_from_map (int imol_map)
 b-factor from map More...
 
SCM map_colour_components (int imol)
 return the colour triple of the imolth map More...
 
PyObject * map_colour_components_py (int imol)
 return the colour triple of the imolth map More...
 
int read_ccp4_map (const std::string &filename, int is_diff_map_flag)
 read a CCP4 map or a CNS map (despite the name).
 
int handle_read_ccp4_map (const std::string &filename, int is_diff_map_flag)
 same function as above - old name for the function. Deleted from the API at some stage
 
int handle_read_emdb_data (const std::string &dir_name)
 
Multi-Residue Torsion
SCM multi_residue_torsion_fit_scm (int imol, SCM residues_specs_scm, int n_trials)
 fit residues More...
 
void multi_residue_torsion_fit (int imol, const std::vector< coot::residue_spec_t > &specs, int n_trials)
 
PyObject * multi_residue_torsion_fit_py (int imol, PyObject *residues_specs_py, int n_trials)
 fit residues More...
 
Merge Fragments

merge fragments

int merge_fragments (int imol)
 each fragment is presumed to be in its own chain.
 
Execute Refmac

execute refmac

void execute_refmac_real (std::string pdb_in_filename, std::string pdb_out_filename, std::string mtz_in_filename, std::string mtz_out_filename, std::string cif_lib_filename, std::string fobs_col_name, std::string sigfobs_col_name, std::string r_free_col_name, short int have_sensible_free_r_flag, short int make_molecules_flag, std::string refmac_count_string, int swap_map_colours_post_refmac_flag, int imol_refmac_map, int diff_map_flag, int phase_combine_flag, std::string phib_string, std::string fom_string, std::string ccp4i_project_dir)
 if swap_map_colours_post_refmac_flag is not 1 thenn imol_refmac_map is ignored.
 
std::string refmac_name (int imol)
 the name for refmac More...
 
Dictionary Functions
std::vector< std::string > dictionary_entries ()
 
void debug_dictionary ()
 
std::string SMILES_for_comp_id (const std::string &comp_id)
 
SCM dictionaries_read ()
 return a list of all the dictionaries read
 
SCM cif_file_for_comp_id_scm (const std::string &comp_id)
 
SCM dictionary_entries_scm ()
 
SCM SMILES_for_comp_id_scm (const std::string &comp_id)
 
PyObject * dictionaries_read_py ()
 
PyObject * cif_file_for_comp_id_py (const std::string &comp_id)
 
PyObject * dictionary_entries_py ()
 
PyObject * SMILES_for_comp_id_py (const std::string &comp_id)
 
Restraints Interface

return the monomer restraints for the given monomer_type, return scheme false on "restraints for monomer not found"

SCM monomer_restraints (const char *monomer_type)
 
*SCM set_monomer_restraints (const char *monomer_type, SCM restraints)
 set the monomer restraints of the given monomer_type More...
 
PyObject * monomer_restraints_py (std::string monomer_type)
 
PyObject * monomer_restraints_for_molecule_py (std::string monomer_type, int imol)
 
PyObject * set_monomer_restraints_py (const char *monomer_type, PyObject *restraints)
 
Atom Information functions

output atom info in a scheme list for use in scripting

in this format (list occ temp-factor element x y z). Return empty list if atom not found.

*SCM atom_info_string_scm (int imol, const char *chain_id, int resno, const char *ins_code, const char *atname, const char *altconf)
 
SCM molecule_to_pdb_string_scm (int imol)
 
std::string resname_from_serial_number (int imol, const char *chain_id, int serial_num)
 return the rename from a residue serial number More...
 
std::string residue_name (int imol, const std::string &chain_id, int resno, const std::string &ins_code)
 return the residue name of the specified residue
 
int serial_number_from_residue_specs (int imol, const std::string &chain_id, int res_no, const std::string &ins_code)
 return the serial number of the specified residue More...
 
SCM residue_info (int imol, const char *chain_id, int resno, const char *ins_code)
 Return a list of atom info for each atom in the specified residue. More...
 
SCM residue_name_scm (int imol, const char *chain_id, int resno, const char *ins_code)
 
SCM chain_fragments_scm (int imol, short int screen_output_also)
 chain fragments
 
int add_molecule (SCM molecule_expression, const char *name)
 generate a molecule from an s-expression More...
 
int clear_and_update_molecule (int molecule_number, SCM molecule_expression)
 update a molecule from a s-expression More...
 
SCM active_residue ()
 return specs of the atom close to screen centre More...
 
SCM closest_atom_simple_scm ()
 return the specs of the closest displayed atom More...
 
SCM closest_atom (int imol)
 return the specs of the closest atom in imolth molecule More...
 
SCM closest_atom_raw_scm ()
 return the specs of the closest atom to the centre of the screen More...
 
SCM residues_near_residue (int imol, SCM residue_in_scm, float radius)
 return residues near residue More...
 
SCM residues_near_residues_scm (int imol, SCM residues_in, float radius)
 return residues near the given residues More...
 
SCM residues_near_position_scm (int imol, SCM pos, float radius)
 residues near residue More...
 
void label_closest_atoms_in_neighbour_residues_scm (int imol, SCM residue_spec_scm, float radius)
 label the closest atoms in the residues that neighbour residue_spec
 
void hydrogenate_region (float radius)
 find the active residue, find the near residues (within radius) create a new molecule, run reduce on that, import hydrogens from the result and apply them to the molecule of the active residue.
 
void add_hydrogens_from_file (int imol, std::string pdb_with_Hs_file_name)
 Add hydrogens to imol from the given pdb file.
 
*PyObject * atom_info_string_py (int imol, const char *chain_id, int resno, const char *ins_code, const char *atname, const char *altconf)
 output atom info in a python list for use in scripting: More...
 
PyObject * molecule_to_pdb_string_py (int imol)
 Return the molecule as a PDB string.
 
PyObject * residue_info_py (int imol, const char *chain_id, int resno, const char *ins_code)
 Return a list of atom info for each atom in the specified residue: More...
 
PyObject * residue_name_py (int imol, const char *chain_id, int resno, const char *ins_code)
 
PyObject * residue_centre_from_spec_py (int imol, PyObject *spec_py)
 
PyObject * chain_fragments_py (int imol, short int screen_output_also)
 
void set_b_factor_residues_py (int imol, PyObject *residue_specs_b_value_tuple_list_py)
 
void set_b_factor_residues_scm (int imol, SCM residue_specs_b_value_tuple_list_scm)
 
Using S-expression molecules
int clear_and_update_molecule_py (int molecule_number, PyObject *molecule_expression)
 
int add_molecule_py (PyObject *molecule_expression, const char *name)
 
PyObject * active_residue_py ()
 Return a list of [imol, chain-id, resno, ins-code, atom-name, alt-conf] for atom that is closest to the screen centre. If there are multiple models with the same coordinates at the screen centre, return the attributes of the atom in the highest number molecule number. More...
 
PyObject * closest_atom_simple_py ()
 return the spec of the closest displayed atom More...
 
PyObject * closest_atom_py (int imol)
 return closest atom in imolth molecule More...
 
PyObject * closest_atom_raw_py ()
 return the specs of the closest atom to the centre of the screen More...
 
PyObject * residues_near_residue_py (int imol, PyObject *residue_in, float radius)
 
PyObject * residues_near_residues_py (int imol, PyObject *residues_in, float radius)
 
PyObject * residues_near_position_py (int imol, PyObject *pos_in, float radius)
 Return residue specs for residues that have atoms that are closer than radius Angstroems to the given position.
 
void label_closest_atoms_in_neighbour_residues_py (int imol, PyObject *residue_spec_py, float radius)
 label the closest atoms in the residues that neighbour residue_spec
 
PyObject * get_bonds_representation (int imol)
 return a Python object for the bonds
 
PyObject * get_dictionary_radii ()
 return a Python object for the radii of the atoms in the dictionary
 
PyObject * get_environment_distances_representation_py (int imol, PyObject *residue_spec_py)
 return a Python object for the representation of bump and hydrogen bonds of
 
PyObject * get_intermediate_atoms_bonds_representation ()
 return a Python object for the intermediate atoms bonds
 
int get_continue_updating_refinement_atoms_state ()
 return the continue-updating-refinement-atoms state
 
status bar string functions
std::string atom_info_as_text_for_statusbar (int atom_index, int imol)
 
std::string atom_info_as_text_for_statusbar (int atom_index, int imol, const std::pair< symm_trans_t, Cell_Translation > &sts)
 
Refinement with specs
SCM all_residues_with_serial_numbers_scm (int imol)
 a utility to return the specs of all the residues, each spec prefixed by the serial number
 
PyObject * all_residues_with_serial_numbers_py (int imol)
 
void regularize_residues (int imol, const std::vector< coot::residue_spec_t > &residues)
 regularize the given residues
 
std::string mtz_file_name (int imol)
 presumes that imol_Refinement_Map has been set
 
SCM refine_zone_with_full_residue_spec_scm (int imol, const char *chain_id, int resno1, const char *inscode_1, int resno2, const char *inscode_2, const char *altconf)
 Refine the given residue range.
 
PyObject * refine_zone_with_full_residue_spec_py (int imol, const char *chain_id, int resno1, const char *inscode_1, int resno2, const char *inscode_2, const char *altconf)
 
void set_show_intermediate_atoms_rota_markup (short int state)
 
void set_show_intermediate_atoms_rama_markup (short int state)
 
void set_cryo_em_refinement (bool mode)
 
bool get_cryo_em_refinement ()
 
SCM accept_moving_atoms_scm ()
 
PyObject * accept_moving_atoms_py ()
 
void register_post_intermediate_atoms_moved_hook (PyObject *function_name)
 
void set_regenerate_bonds_needs_make_bonds_type_checked (bool state)
 
bool get_regenerate_bonds_needs_make_bonds_type_checked_state ()
 
Water Chain Functions
SCM water_chain_from_shelx_ins_scm (int imol)
 return the chain id of the water chain from a shelx molecule. Raw interface More...
 
SCM water_chain_scm (int imol)
 return the chain id of the water chain. Raw interface
 
PyObject * water_chain_from_shelx_ins_py (int imol)
 
PyObject * water_chain_py (int imol)
 return the chain id of the water chain. Raw interface
 
Glyco Tools
void print_glyco_tree (int imol, const std::string &chain_id, int resno, const std::string &ins_code)
 print the glycosylation tree that contains the specified residue
 
Variance Map

Make a variance map, based on the grid of the first map.

Returns
the molecule number of the new map. Return -1 if unable to make a variance map.
int make_variance_map (std::vector< int > map_molecule_number_vec)
 
int make_variance_map_scm (SCM map_molecule_number_list)
 
int make_variance_map_py (PyObject *map_molecule_number_list)
 
Spin Search Functions
void spin_search_by_atom_vectors (int imol_map, int imol, const std::string &chain_id, int resno, const std::string &ins_code, const std::pair< std::string, std::string > &direction_atoms_list, const std::vector< std::string > &moving_atoms_list)
 
void spin_search (int imol_map, int imol, const char *chain_id, int resno, const char *ins_code, SCM direction_atoms_list, SCM moving_atoms_list)
 for the given residue, spin the atoms in moving_atom_list around the bond defined by direction_atoms_list looking for the best fit to density of imol_map map of the first atom in moving_atom_list. Works (only) with atoms in altconf ""
 
void spin_N_scm (int imol, SCM residue_spec_scm, float angle)
 Spin N and CB (and the rest of the side chain if extant) More...
 
SCM CG_spin_search_scm (int imol_model, int imol_map)
 Spin search the density based on possible positions of CG of a side-chain. More...
 
void spin_search_py (int imol_map, int imol, const char *chain_id, int resno, const char *ins_code, PyObject *direction_atoms_list, PyObject *moving_atoms_list)
 for the given residue, spin the atoms in moving_atom_list... More...
 
void spin_N_py (int imol, PyObject *residue_spec, float angle)
 Spin N and CB (and the rest of the side chain if extant) More...
 
PyObject * CG_spin_search_py (int imol_model, int imol_map)
 Spin search the density based on possible positions of CG of a side-chain.
 
Extra Ligand Functions

make conformers of the ligand search molecules, each in its own molecule.

Don't search the density.

! !

Return a list of new molecule numbers

SCM ligand_search_make_conformers_scm ()
 
PyObject * ligand_search_make_conformers_py ()
 
std::vector< int > ligand_search_make_conformers_internal ()
 
Dock Sidechains

cootaneer (i.e. assign sidechains onto mainchain model)

atom_in_fragment_atom_spec is any atom spec in the fragment that should be assigned with sidechains.

Returns
the success status (0 is fail).
int cootaneer (int imol_map, int imol_model, SCM atom_in_fragment_atom_spec)
 
int cootaneer_py (int imol_map, int imol_model, PyObject *atom_in_fragment_atom_spec)
 
Sequence from Map

Use the map to estimate the sequence - you will need a decent map

Returns
the guessed sequence (empty is fail).
std::string sequence_from_map (int imol, const std::string &chain_id, int resno_start, int resno_end, int imol_map)
 
void apply_sequence_to_fragment (int imol, const std::string &chain_id, int resno_start, int resno_end, int imol_map, const std::string &file_name_for_sequences)
 
void assign_sequence_to_active_fragment ()
 
Rotamer Scoring
std::vector
< coot::named_rotamer_score > 
score_rotamers (int imol, const char *chain_id, int res_no, const char *ins_code, const char *alt_conf, int imol_map, int clash_flag, float lowest_probability)
 
SCM score_rotamers_scm (int imol, const char *chain_id, int res_no, const char *ins_code, const char *alt_conf, int imol_map, int clash_flag, float lowest_probability)
 return the scores of the rotamers for this residue.
 
PyObject * score_rotamers_py (int imol, const char *chain_id, int res_no, const char *ins_code, const char *alt_conf, int imol_map, int clash_flag, float lowest_probability)
 
protein-db
std::pair< std::pair< int, int >
, std::vector< int > > 
protein_db_loops (int imol_coords, const std::vector< coot::residue_spec_t > &residue_specs, int imol_map, int nfrags, bool preserve_residue_names)
 Cowtan's protein_db loops.
 
std::string protein_db_loop_specs_to_atom_selection_string (const std::vector< coot::residue_spec_t > &specs)
 
SCM protein_db_loops_scm (int imol_coords, SCM residues_specs, int imol_map, int nfrags, bool preserve_residue_names)
 
PyObject * protein_db_loops_py (int imol_coords, PyObject *residues_specs, int imol_map, int nfrags, bool preserve_residue_names)
 
Coot's Hole implementation
void hole (int imol, float start_x, float start_y, float start_z, float end_x, float end_y, float end_z, float colour_map_multiplier, float colour_map_offset, int n_runs, bool show_probe_radius_graph_flag, std::string export_surface_dots_file_name)
 starting piont and end point, colour map multiplier and shall the probe radius graph be shown (dummy value currently). More...
 
void probe_radius_graph_close_callback (GtkWidget *button, GtkWidget *dialog)
 
void show_hole_probe_radius_graph (const std::vector< std::pair< clipper::Coord_orth, double > > &hole_path, double path_length)
 
void show_hole_probe_radius_graph_basic (const std::vector< std::pair< clipper::Coord_orth, double > > &hole_path, double path_length)
 
void show_hole_probe_radius_graph_goocanvas (const std::vector< std::pair< clipper::Coord_orth, double > > &hole_path, double path_length)
 
void make_link (int imol, coot::atom_spec_t &spec_1, coot::atom_spec_t &spec_2, const std::string &link_name, float length)
 make a link between the specified atoms
 
void make_link_scm (int imol, SCM spec_1, SCM spec_2, const std::string &link_name, float length)
 
SCM link_info_scm (int imol)
 
void make_link_py (int imol, PyObject *spec_1, PyObject *spec_2, const std::string &link_name, float length)
 
PyObject * link_info_py (int imol)
 
Drag and Drop Functions
int handle_drag_and_drop_string (const std::string &uri)
 handle the string that get when a file or URL is dropped.
 
Map Contouring Functions
PyObject * map_contours (int imol, float contour_level)
 return two lists: a list of vertices and a list of indices for connection
 
Map to Model Correlation
void set_map_correlation_atom_radius (float r)
 The atom radius is not passed as a parameter to correlation.
 
SCM map_to_model_correlation_scm (int imol, SCM residue_specs, SCM neighb_residue_specs, unsigned short int atom_mask_mode, int imol_map)
 atom-mask-mode is as follows:
 
SCM map_to_model_correlation_stats_scm (int imol, SCM residue_specs, SCM neighb_residue_specs, unsigned short int atom_mask_mode, int imol_map)
 
PyObject * map_to_model_correlation_py (int imol, PyObject *residue_specs, PyObject *neighb_residue_specs, unsigned short int atom_mask_mode, int imol_map)
 
PyObject * map_to_model_correlation_stats_py (int imol, PyObject *residue_specs, PyObject *neighb_residue_specs, unsigned short int atom_mask_mode, int imol_map)
 
PyObject * map_to_model_correlation_stats_per_residue_range_py (int imol, const std::string &chain_id, int imol_map, unsigned int n_residue_per_residue_range, short int exclude_NOC_flag)
 
float map_to_model_correlation (int imol, const std::vector< coot::residue_spec_t > &residue_specs, const std::vector< coot::residue_spec_t > &neigh_residue_specs, unsigned short int atom_mask_mode, int imol_map)
 atom-mask-mode is as follows:
 
coot::util::density_correlation_stats_info_t map_to_model_correlation_stats (int imol, const std::vector< coot::residue_spec_t > &residue_specs, const std::vector< coot::residue_spec_t > &neigh_residue_specs, unsigned short int atom_mask_mode, int imol_map)
 map to model density correlation stats More...
 
std::vector< std::pair
< coot::residue_spec_t, float > > 
map_to_model_correlation_per_residue (int imol, const std::vector< coot::residue_spec_t > &specs, unsigned short int atom_mask_mode, int imol_map)
 map to model density correlation, reported per residue More...
 
std::map< coot::residue_spec_t,
coot::util::density_stats_info_t > 
map_to_model_correlation_stats_per_residue (int imol, const std::vector< coot::residue_spec_t > &residue_specs, unsigned short int atom_mask_mode, float atom_radius_for_masking, int imol_map)
 map to model density statistics, reported per residue
 
std::pair< std::map
< coot::residue_spec_t,
coot::util::density_correlation_stats_info_t >
, std::map
< coot::residue_spec_t,
coot::util::density_correlation_stats_info_t > > 
map_to_model_correlation_stats_per_residue_range (int imol, const std::string &chain_id, int imol_map, unsigned int n_residue_per_residue_range, short int exclude_NOC_flag)
 map to model density statistics, reported per residue, the middle residue of a range of residues More...
 
SCM map_to_model_correlation_per_residue_scm (int imol, SCM residue_specs, unsigned short int atom_mask_mode, int imol_map)
 map to model correlation
 
SCM map_to_model_correlation_stats_per_residue_scm (int imol, SCM residue_specs_scm, unsigned short int atom_mask_mode, float atom_radius_for_masking, int imol_map)
 map to model stats
 
SCM map_to_model_correlation_stats_per_residue_range_scm (int imol, const std::string &chain_id, int imol_map, unsigned int n_residue_per_residue_range, short int exclude_NOC_flag)
 
SCM qq_plot_map_and_model_scm (int imol, SCM residue_specs_scm, SCM neigh_residue_specs_scm, unsigned short int atom_mask_mode, int imol_map)
 QQ plot of the model density correlation, reported per residue. More...
 
PyObject * map_to_model_correlation_per_residue_py (int imol, PyObject *residue_specs, unsigned short int atom_mask_mode, int imol_map)
 
PyObject * qq_plot_map_and_model_py (int imol, PyObject *residue_specs_py, PyObject *neigh_residue_specs_py, unsigned short int atom_mask_mode, int imol_map)
 
float density_score_residue_scm (int imol, SCM residue_spec, int imol_map)
 
float density_score_residue_py (int imol, PyObject *residue_spec, int imol_map)
 
float density_score_residue (int imol, const char *chain_id, int res_no, const char *ins_code, int imol_map)
 simple density score for given residue (over-ridden by scripting function)
 
SCM map_mean_scm (int imol)
 return sigma for the given map. Return scheme False if not a valid map molecule number.
 
SCM map_sigma_scm (int imol)
 
SCM map_statistics_scm (int imol)
 return either scheme false on non-a-map or list (mean, standard-deviation, skew, kurtosis)
 
PyObject * map_mean_py (int imol)
 return sigma for the given map. Return Python False if not a valid map molecule number.
 
PyObject * map_sigma_py (int imol)
 
PyObject * map_statistics_py (int imol)
 
Get Sequence
std::string get_sequence_as_fasta_for_chain (int imol, const std::string &chain_id)
 get the sequence for chain_id in imol
 
void write_sequence (int imol, const std::string &file_name)
 write the sequence for imol as fasta
 

Detailed Description

Coot Scripting Interface - General (C++ functions)

Function Documentation

SCM active_residue ( )

return specs of the atom close to screen centre

Return a list of (list imol chain-id resno ins-code atom-name alt-conf) for atom that is closest to the screen centre in any displayed molecule. If there are multiple models with the same coordinates at the screen centre, return the attributes of the atom in the highest number molecule number.

return scheme false if no active residue

PyObject* active_residue_py ( )

Return a list of [imol, chain-id, resno, ins-code, atom-name, alt-conf] for atom that is closest to the screen centre. If there are multiple models with the same coordinates at the screen centre, return the attributes of the atom in the highest number molecule number.

return False if no active residue

int add_molecule ( SCM  molecule_expression,
const char *  name 
)

generate a molecule from an s-expression

return a molecule number, -1 on error

std::pair<int, std::string> align_to_closest_chain ( std::string  target_seq,
float  match_fraction 
)

align sequence to closest chain (compare across all chains in all molecules).

Typically match_fraction is 0.95 or so.

Return the molecule number and chain id if successful, return -1 as the molecule number if not.

* PyObject* atom_info_string_py ( int  imol,
const char *  chain_id,
int  resno,
const char *  ins_code,
const char *  atname,
const char *  altconf 
)

output atom info in a python list for use in scripting:

in this format [occ, temp_factor, element, x, y, z]. Return empty list if atom not found.

std::vector<int> auto_read_make_and_draw_maps ( const char *  filename)

read MTZ file filename and from it try to make maps

Useful for reading the output of refmac. The default labels (FWT/PHWT and DELFWT/PHDELFWT) can be changed using ...[something]

Returns
a list of molecule numbers for the new maps
float b_factor_from_map ( int  imol_map)

b-factor from map

calculate structure factors and use the amplitudes to estimate the B-factor of the data using a wilson plot using a low resolution limit of 4.5A.

Returns
-1 when given a bad map or there were no data beyond 4.5A
SCM CG_spin_search_scm ( int  imol_model,
int  imol_map 
)

Spin search the density based on possible positions of CG of a side-chain.

c.f. EM-Ringer

int clashes_with_symmetry ( int  imol,
const char *  chain_id,
int  res_no,
const char *  ins_code,
float  clash_dist 
)

return 1 if this residue clashes with the symmetry-related atoms of the same molecule.

0 means that it did not clash, -1 means that the residue or molecule could not be found or that there was no cell and symmetry.

int clear_and_update_molecule ( int  molecule_number,
SCM  molecule_expression 
)

update a molecule from a s-expression

And going the other way, given an s-expression, update molecule_number by the given molecule. Clear what's currently there first though.

SCM closest_atom ( int  imol)

return the specs of the closest atom in imolth molecule

Return a list of (list imol chain-id resno ins-code atom-name alt-conf (list x y z)) for atom that is closest to the screen centre in the given molecule (unlike active-residue, no account is taken of the displayed state of the molecule). If there is no atom, or if imol is not a valid model molecule, return scheme false.

PyObject* closest_atom_py ( int  imol)

return closest atom in imolth molecule

Return a list of [imol, chain-id, resno, ins-code, atom-name, alt-conf, [x, y, z]] for atom that is closest to the screen centre in the given molecule (unlike active-residue, no account is taken of the displayed state of the molecule). If there is no atom, or if imol is not a valid model molecule, return False.

PyObject* closest_atom_raw_py ( )

return the specs of the closest atom to the centre of the screen

Return a list of (list imol chain-id resno ins-code atom-name alt-conf (list x y z)) for atom that is closest to the screen for displayed molecules. If there is no atom, return scheme false. Don't choose the CA of the residue if there is a CA in the residue of the closest atom

SCM closest_atom_raw_scm ( )

return the specs of the closest atom to the centre of the screen

Return a list of (list imol chain-id resno ins-code atom-name alt-conf (list x y z)) for atom that is closest to the screen for displayed molecules. If there is no atom, return scheme false. Don't choose the CA of the residue if there is a CA in the residue of the closest atom. 201602015-PE: I add this now, but I have a feeling that I've done this before.

PyObject* closest_atom_simple_py ( )

return the spec of the closest displayed atom

Return a list of [imol, chain-id, resno, ins-code, atom-name, alt-conf, [x, y, z]] for atom that is closest to the screen centre in the given molecule (unlike active-residue, potential CA substition is not performed). If there is no atom, or if imol is not a valid model molecule, return False.

SCM closest_atom_simple_scm ( )

return the specs of the closest displayed atom

Return a list of (list imol chain-id resno ins-code atom-name alt-conf (list x y z)) for atom that is closest to the screen centre in the given molecule (unlike active-residue, potential CA substition is not performed). If there is no atom, or if imol is not a valid model molecule, return scheme false.

int flip_hand ( int  imol_map)

Flip the hand of the map.

in case it was accidentally generated on the wrong one.

Returns
the molecule number of the flipped map.
void get_coords_for_accession_code ( const std::string &  code)

if possible, read in the new coords getting coords via web.

(no return value because get-url-str does not return one).

* SCM get_symmetry ( int  imol)

return the symmetry of the imolth molecule

Return as a list of strings the symmetry operators of the given molecule. If imol is a not a valid molecule, return an empty list.

void go_to_map_molecule_centre ( int  imol_map)

Go to the centre of the molecule - for Cryo-EM Molecules.

and recontour at a sensible value.

void hole ( int  imol,
float  start_x,
float  start_y,
float  start_z,
float  end_x,
float  end_y,
float  end_z,
float  colour_map_multiplier,
float  colour_map_offset,
int  n_runs,
bool  show_probe_radius_graph_flag,
std::string  export_surface_dots_file_name 
)

starting piont and end point, colour map multiplier and shall the probe radius graph be shown (dummy value currently).

if export_dots_file_name string length is zero, then don't try to export the surface dots.

SCM map_colour_components ( int  imol)

return the colour triple of the imolth map

(e.g.: (list 0.4 0.6 0.8). If invalid imol return scheme false.

PyObject* map_colour_components_py ( int  imol)

return the colour triple of the imolth map

e.g.: [0.4, 0.6, 0.8]. If invalid imol return Py_False.

std::vector<std::pair<coot::residue_spec_t,float> > map_to_model_correlation_per_residue ( int  imol,
const std::vector< coot::residue_spec_t > &  specs,
unsigned short int  atom_mask_mode,
int  imol_map 
)

map to model density correlation, reported per residue

atom-mask-mode is as follows:

coot::util::density_correlation_stats_info_t map_to_model_correlation_stats ( int  imol,
const std::vector< coot::residue_spec_t > &  residue_specs,
const std::vector< coot::residue_spec_t > &  neigh_residue_specs,
unsigned short int  atom_mask_mode,
int  imol_map 
)

map to model density correlation stats

atom-mask-mode is as follows:

std::pair<std::map<coot::residue_spec_t, coot::util::density_correlation_stats_info_t>, std::map<coot::residue_spec_t, coot::util::density_correlation_stats_info_t> > map_to_model_correlation_stats_per_residue_range ( int  imol,
const std::string &  chain_id,
int  imol_map,
unsigned int  n_residue_per_residue_range,
short int  exclude_NOC_flag 
)

map to model density statistics, reported per residue, the middle residue of a range of residues

Returns
the all-atom stats first and side chains stats second
PyObject* multi_residue_torsion_fit_py ( int  imol,
PyObject *  residues_specs_py,
int  n_trials 
)

fit residues

(note: fit to the current-refinement map)

SCM multi_residue_torsion_fit_scm ( int  imol,
SCM  residues_specs_scm,
int  n_trials 
)

fit residues

(note: fit to the current-refinement map)

void multi_sharpen_blur_map_py ( int  imol_map,
PyObject *  b_factors_list 
)

make many sharpened or blurred maps

blurred maps are generated by using a positive value of b_factor.

void multi_sharpen_blur_map_scm ( int  imol_map,
SCM  b_factors_list 
)

make many sharpened or blurred maps

blurred maps are generated by using a positive value of b_factor.

void prodrg_import_function ( std::string  file_name,
std::string  comp_id 
)

import given mdl file into prodrg or other 3d generation program

the function passed to lbg, so that it calls it when a new prodrg-in.mdl file has been made. We no longer have a timeout function waiting for prodrg-in.mdl to be updated/written.

SCM qq_plot_map_and_model_scm ( int  imol,
SCM  residue_specs_scm,
SCM  neigh_residue_specs_scm,
unsigned short int  atom_mask_mode,
int  imol_map 
)

QQ plot of the model density correlation, reported per residue.

atom-mask-mode is as follows:

std::string refmac_name ( int  imol)

the name for refmac

Returns
a stub name used in the construction of filename for refmac output
SCM residue_info ( int  imol,
const char *  chain_id,
int  resno,
const char *  ins_code 
)

Return a list of atom info for each atom in the specified residue.

output is like this: (list (list (list atom-name alt-conf) (list occ temp-fact element) (list x y z)))

occ can be a single number or a list of seven numbers of which the first is the isotropic B.

PyObject* residue_info_py ( int  imol,
const char *  chain_id,
int  resno,
const char *  ins_code 
)

Return a list of atom info for each atom in the specified residue:

output is like this: [ [[atom-name,alt-conf] [occ,temp_fact,element] [x,y,z]]]

SCM residues_near_position_scm ( int  imol,
SCM  pos,
float  radius 
)

residues near residue

Returns
residues within radius of pos (x,y,z) position

Return a list of pairs of (imol, residue_spec). pos is a list of 3 numbers. (get imol from active-atom)

SCM residues_near_residue ( int  imol,
SCM  residue_in_scm,
float  radius 
)

return residues near residue

Return residue specs for residues that have atoms that are closer than radius Angstroems to any atom in the residue specified by res_in.

SCM residues_near_residues_scm ( int  imol,
SCM  residues_in,
float  radius 
)

return residues near the given residues

Return residue specs for residues that have atoms that are closer than radius Angstroems to any atom in the residue specified by res_in.

std::string resname_from_serial_number ( int  imol,
const char *  chain_id,
int  serial_num 
)

return the rename from a residue serial number

Returns
blank ("") on failure.
void sbase_import_function ( std::string  comp_id)

import molecule from CCP4 SRS (or SBase, as it used to be called).

the function passed to lbg, so that it calls it when a new SBase comp_id is required. We no longer have a timeout function waiting for prodrg-in.mdl to be updated/written.

int serial_number_from_residue_specs ( int  imol,
const std::string &  chain_id,
int  res_no,
const std::string &  ins_code 
)

return the serial number of the specified residue

Returns
-1 on failure to find the residue
* SCM set_monomer_restraints ( const char *  monomer_type,
SCM  restraints 
)

set the monomer restraints of the given monomer_type

Returns
scheme false or true for success or failure to set the restrains for monomer_type
void set_radial_map_colouring_saturation ( int  imol,
float  saturation 
)

radial map colouring saturation

saturation is a number between 0 and 1, typically 0.5

int sharpen_blur_map ( int  imol_map,
float  b_factor 
)

make a sharpened or blurred map

blurred maps are generated by using a positive value of b_factor.

Returns
the index of the map created by applying a b-factor to the given map. Return -1 on failure.
int sharpen_blur_map_with_resampling ( int  imol_map,
float  b_factor,
float  resample_factor 
)

make a sharpened or blurred map with resampling

resampling factor might typically be 1.3

blurred maps are generated by using a positive value of b_factor.

Returns
the index of the map created by applying a b-factor to the given map. Return -1 on failure.
void spin_N_py ( int  imol,
PyObject *  residue_spec,
float  angle 
)

Spin N and CB (and the rest of the side chain if extant)

Sometime on N-terminal addition, then N ends up pointing the wrong way. The allows us to (more or less) interchange the positions of the CB and the N. angle is in degrees.

void spin_N_scm ( int  imol,
SCM  residue_spec_scm,
float  angle 
)

Spin N and CB (and the rest of the side chain if extant)

Sometime on N-terminal addition, then N ends up pointing the wrong way. The allows us to (more or less) interchange the positions of the CB and the N. angle is in degrees.

void spin_search_py ( int  imol_map,
int  imol,
const char *  chain_id,
int  resno,
const char *  ins_code,
PyObject *  direction_atoms_list,
PyObject *  moving_atoms_list 
)

for the given residue, spin the atoms in moving_atom_list...

around the bond defined by direction_atoms_list looking for the best fit to density of imom_map map of the first atom in moving_atom_list. Works (only) with atoms in altconf ""

SCM water_chain_from_shelx_ins_scm ( int  imol)

return the chain id of the water chain from a shelx molecule. Raw interface

Returns
scheme false if no chain or bad imol