From bef6c7548e536cbd64d3d3d422f1c29bf9ba7e28 Mon Sep 17 00:00:00 2001
From: Egor Sobolev <egor.sobolev@xfel.eu>
Date: Fri, 28 Jul 2023 13:23:57 +0200
Subject: [PATCH] Add cell volume calculation and reading of unit cell file
---
src/geomtools/sfx/lattice.py | 56 ++++++++++++++++++++++++++++++++----
1 file changed, 50 insertions(+), 6 deletions(-)
diff --git a/src/geomtools/sfx/lattice.py b/src/geomtools/sfx/lattice.py
index 92c8c26..20059f8 100644
--- a/src/geomtools/sfx/lattice.py
+++ b/src/geomtools/sfx/lattice.py
@@ -2,7 +2,6 @@
import numpy as np
-
# Electron charge in C
ELECTRON_CHARGE = 1.6021773e-19
# Planck's constant (Js)
@@ -43,13 +42,58 @@ def get_min_bragg_dist(r, clen, lmd, cell):
return min(invq / a for a in cell)
+def read_crystfel_cell(stream):
+ state = 0
+ cell = {}
+ for l in stream:
+ ln, _, _ = l.partition(';')
+ ln = ln.strip()
+ if state == 0:
+ if ln and ln.startswith("CrystFEL unit cell file"):
+ state = 1
+ elif ln:
+ raise ValueError("This is not CrystFEL unit cell file")
+ elif state == 1:
+ if ln:
+ name, _, val = ln.partition('=')
+ name = name.strip()
+ val = val.strip()
+ if name in ['a', 'b', 'c']:
+ num, _, unit = val.partition(' ')
+ num = float(num)
+ if unit == 'nm':
+ num *= 10
+ elif name in ['al', 'be', 'ga']:
+ num, _, unit = val.partition(' ')
+ num = float(num)
+ if unit == 'rad':
+ num = num / np.pi * 180
+ else:
+ num = val
+
+ cell[name] = num
+
+ return type('UnitCell', (), cell)
+
+
+def cell_volume(a, b, c, alpha, beta, gamma):
+ al = alpha * np.pi / 180.
+ be = beta * np.pi / 180.
+ ga = gamma * np.pi / 180.
+ abc = a * b * c
+ csa, csb, csc = np.cos(al), np.cos(be), np.cos(ga)
+ V2 = abc * abc * (
+ 1. - csa * csa - csb * csb - csc * csc + 2. * csa * csb * csc)
+ return np.sqrt(V2)
+
+
def spacing(h, k, l, a, b, c, alpha, beta, gamma):
- alpha *= np.pi / 180.
- beta *= np.pi / 180.
- gamma *= np.pi / 180.
+ al = alpha * np.pi / 180.
+ be = beta * np.pi / 180.
+ ga = gamma * np.pi / 180.
a2, b2, c2, abc = a * a, b * b, c * c, a * b * c
- sna, snb, snc = np.sin(alpha), np.sin(beta), np.sin(gamma)
- csa, csb, csc = np.cos(alpha), np.cos(beta), np.cos(gamma)
+ sna, snb, snc = np.sin(al), np.sin(be), np.sin(ga)
+ csa, csb, csc = np.cos(al), np.cos(be), np.cos(ga)
V2 = abc * abc * (
1. - csa * csa - csb * csb - csc * csc + 2. * csa * csb * csc)
S11 = b2 * c2 * sna * sna * h * h
--
GitLab