Skip to content

Commit 40ce9b6

Browse files
committed
Tidied Files NO_JIRA
1 parent 133920b commit 40ce9b6

File tree

3 files changed

+177
-0
lines changed

3 files changed

+177
-0
lines changed

assets/save_script.gif

17.7 MB
Loading
Lines changed: 25 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,25 @@
1+
# Analyse Disorder and Geometry
2+
3+
## Summary
4+
5+
This script demonstrates how to analyse crystallographic disorder and evaluate molecular geometry using the Cambridge Structural Database (CSD) Python API. It inspects disorder groups in a crystal structure, activates specific disorder configurations, identifies unusual torsion angles, and visualises torsion histograms using Matplotlib.
6+
7+
## Requirements
8+
9+
CSD Python API (v 3.5.0 or later)
10+
11+
Matplotlib (v 3.9.4 or later)
12+
13+
## Licensing Requirements
14+
15+
CSD-Core
16+
17+
## Instructions on Running
18+
19+
In a command prompt run `python analyse_disorder.py`.
20+
21+
## Author
22+
23+
Andrew J. Peel (2025-09-29)
24+
25+
> For feedback or to report any issues please contact [support@ccdc.cam.ac.uk](mailto:support@ccdc.cam.ac.uk)
Lines changed: 152 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,152 @@
1+
#! /usr/bin/env python
2+
###############################################################################
3+
#
4+
# This script can be used for any purpose without limitation subject to the
5+
# conditions at http://www.ccdc.cam.ac.uk/Community/Pages/Licences/v2.aspx
6+
#
7+
# This permission notice and the following statement of attribution must be
8+
# included in all copies or substantial portions of this script.
9+
#
10+
# 2025-09-29: created by Andrew J. Peel, Cambridge Crystallographic Data Centre
11+
#
12+
###############################################################################
13+
14+
from ccdc.io import EntryReader
15+
from ccdc.conformer import GeometryAnalyser
16+
import matplotlib.pyplot as plt
17+
18+
# Read CSD entry
19+
entry_reader = EntryReader('CSD')
20+
purzio_entry = entry_reader.entry('PURZIO')
21+
print(f'Refcode: {purzio_entry.identifier}')
22+
23+
# Get crystal
24+
purzio_crystal = purzio_entry.crystal
25+
26+
# Examine disorder groups
27+
print(f'Crystal has disorder: {purzio_crystal.has_disorder}')
28+
print(f'Number of disorder assemblies: '
29+
f'{len(purzio_crystal.disorder.assemblies)}')
30+
for i, assembly in enumerate(purzio_crystal.disorder.assemblies, start=1):
31+
print(f'Assembly {i} has {len(assembly.groups)} disorder groups.')
32+
33+
# Activate disorder assembly, group 1 of 2
34+
first_assembly = purzio_crystal.disorder.assemblies[0]
35+
first_group = first_assembly.groups[0]
36+
first_group.activate()
37+
38+
print(f'Activated assembly_id: {first_assembly.id}')
39+
print(f'Activated group_id: {first_assembly.active.id}')
40+
print(f'Occupancy of disorder group: {first_group.occupancy:.3f}')
41+
42+
# Examine the molecules
43+
mol = purzio_crystal.molecule
44+
print(f'Number of molecules in crystal: '
45+
f'{len( mol.components)}')
46+
47+
# Select component by index
48+
comp_index = 0
49+
50+
# Examine the first molecule
51+
comp = mol.components[comp_index]
52+
print(f'Analysing molecule {comp_index + 1} of '
53+
f'{len(mol.components)}')
54+
55+
# Check if this molecule is disordered
56+
if any([atom.occupancy < 1 for atom in comp.atoms]):
57+
print('This molecule is disordered.')
58+
else:
59+
print('This molecule is not disordered.')
60+
61+
# Identify disordered atoms
62+
print([(atom.label, atom.occupancy) for atom in comp.atoms
63+
if atom.occupancy < 1])
64+
65+
# Create an instance of geometry analyser
66+
engine = GeometryAnalyser()
67+
68+
# Adjust settings - exclude organometallics and powder structures
69+
engine.settings.organometallic_filter = 'organics_only'
70+
engine.settings.powder_filter = True
71+
72+
# Analyse geometry
73+
geom_analysed_mol = engine.analyse_molecule(comp)
74+
75+
# For example, look at torsion measurements
76+
for torsion in geom_analysed_mol.analysed_torsions:
77+
if torsion.unusual and torsion.enough_hits:
78+
print(f'Unusual torsion: {",".join(torsion.atom_labels)} '
79+
f'Torsion: {torsion.value:.3f} °')
80+
81+
# Accessing histogram data
82+
# Take example of an unusual torsion
83+
unusual_tors = [t for t in geom_analysed_mol.analysed_torsions
84+
if t.unusual and t.enough_hits]
85+
86+
# Retrieves the first occurance of specified torsion
87+
tor_query = next(t for t in unusual_tors if
88+
t.atom_labels == ['C39', 'C28', 'C29', 'C30'])
89+
90+
counts = tor_query.histogram()
91+
92+
# Create histogram
93+
n_bins = len(counts)
94+
bin_min, bin_max = 0, 180
95+
bin_width = (bin_max - bin_min) / n_bins
96+
bin_edges = [bin_min + i * bin_width for i in range(n_bins + 1)]
97+
98+
plt.figure(figsize=(8, 5))
99+
plt.bar(
100+
bin_edges[:-1],
101+
counts,
102+
width=bin_width,
103+
edgecolor="black",
104+
align="edge",
105+
color="limegreen"
106+
)
107+
108+
# Red line for torsion value
109+
torsion_val = tor_query.value
110+
plt.axvline(
111+
torsion_val,
112+
color="red",
113+
linestyle="--",
114+
linewidth=2
115+
)
116+
117+
# Add text label above the histogram, slightly shifted
118+
y_max = max(counts)
119+
plt.text(
120+
torsion_val + 2, # shift slightly right
121+
y_max * 0.9, # place near the top
122+
f"{torsion_val:.2f}°",
123+
color="red",
124+
fontsize=10,
125+
rotation=90,
126+
va="bottom"
127+
)
128+
129+
plt.xlabel("Torsion Angle (degrees)")
130+
plt.ylabel("Count of Observation Hits")
131+
plt.title(f"Histogram for torsion {tor_query.atom_labels}")
132+
plt.xlim(bin_min, bin_max)
133+
134+
plt.show()
135+
136+
print('\n----------------------------------\n')
137+
print('Now using disorder combinations\n')
138+
139+
# Interate over disorder combinations
140+
for comb_id, combination in enumerate(purzio_crystal.disorder.combinations):
141+
for comp_id, component in enumerate(purzio_crystal.molecule.components):
142+
print(f'Combination_id {comb_id}, Component_id {comp_id}')
143+
analysed_mol = engine.analyse_molecule(component)
144+
unusual_tors = [t for t in analysed_mol.analysed_torsions
145+
if t.enough_hits and t.unusual]
146+
if len(unusual_tors) > 0:
147+
for t in unusual_tors:
148+
print(f'Unusual torsion: {",".join(t.atom_labels)}'
149+
f' Torsion: {t.value:.3f} °')
150+
else:
151+
print('No unusual torsions with enough hits.')
152+
print()

0 commit comments

Comments
 (0)