diff --git a/MDANSE/Src/MDANSE/Framework/Configurators/CorrelationFramesConfigurator.py b/MDANSE/Src/MDANSE/Framework/Configurators/CorrelationFramesConfigurator.py
index bc341933d1..294038f488 100644
--- a/MDANSE/Src/MDANSE/Framework/Configurators/CorrelationFramesConfigurator.py
+++ b/MDANSE/Src/MDANSE/Framework/Configurators/CorrelationFramesConfigurator.py
@@ -19,15 +19,21 @@
class CorrelationFramesConfigurator(FramesConfigurator):
+ """Parses the input of trajectory frames.
+
+ Configures the time frame range to be used in the calculations
+ together with a movable window used for correlations.
+ """
+
def configure(self, value: tuple[int, int, int, int]):
- """Configure the correlation and set the number of correlation
- frames to use.
+ """Set the number of correlation frames to use.
Parameters
----------
value : tuple[int, int, int, int]
The frames setting plus the number of frames used for the
correlations.
+
"""
trajConfig = self._configurable[self._dependencies["trajectory"]]
n_steps = trajConfig["length"]
@@ -38,6 +44,7 @@ def configure(self, value: tuple[int, int, int, int]):
first, last, step, c_frames = value
super().configure((first, last, step))
+ self._original_input = value
if c_frames > self["n_frames"]:
self.error_status = (
diff --git a/MDANSE/Src/MDANSE/Framework/Configurators/QVectorsConfigurator.py b/MDANSE/Src/MDANSE/Framework/Configurators/QVectorsConfigurator.py
index a9fde6687e..e4b5d1003a 100644
--- a/MDANSE/Src/MDANSE/Framework/Configurators/QVectorsConfigurator.py
+++ b/MDANSE/Src/MDANSE/Framework/Configurators/QVectorsConfigurator.py
@@ -14,28 +14,28 @@
# along with this program. If not, see .
#
+from typing import Any
from MDANSE.Framework.Configurators.IConfigurator import IConfigurator
from MDANSE.Framework.QVectors.IQVectors import IQVectors
class QVectorsConfigurator(IConfigurator):
- """
- This Configurator allows to set reciprocal vectors.
-
- Reciprocal vectors are used in MDANSE for the analysis related to scattering experiments such as dynamic coherent structure
- or elastic incoherent structure factor analysis. In MDANSE, properties that depends on Q vectors are always scalar regarding
- Q vectors in the sense that the values of these properties will be computed for a given norm of Q vectors and not for a given Q vectors.
- Hence, the Q vectors generator supported by MDANSE always generates Q vectors on Q-shells, each shell containing a set of Q vectors whose
+ """Creates and configures a q-vector generator.
+
+ Reciprocal vectors are used in MDANSE for analysis related to
+ scattering experiments, such as dynamic coherent structure
+ or elastic incoherent structure factor analysis. In MDANSE, properties
+ that depend on Q vectors are always scalar regarding Q vectors
+ in the sense that the values of these properties will be computed
+ for a given norm of Q vectors and not for a given Q vector.
+ Hence, the Q vectors generator supported by MDANSE always generates
+ Q vectors on Q-shells, each shell containing a set of Q vectors whose
norm match the Q shell value within a given tolerance.
- Depending on the generator selected, Q vectors can be generated isotropically or anistropically, on a lattice or randomly.
-
- Q vectors can be saved to a user definition and, as such, can be further reused in another MDANSE session.
+ Depending on the generator selected, Q vectors can be generated
+ isotropically or anistropically, on a lattice or randomly.
- To define a new Q vectors generator, you must inherit from MDANSE.Framework.QVectors.QVectors.QVector interface.
-
- :note: this configurator depends on 'trajectory' configurator to be configured.
"""
_default = (
@@ -43,27 +43,26 @@ class QVectorsConfigurator(IConfigurator):
{"shells": (0.1, 5, 0.1), "width": 0.1, "n_vectors": 50, "seed": 0},
)
- def configure(self, value):
- """
- Configure a Q vectors generator.
+ def configure(self, value: tuple[str, dict[str, Any]]):
+ """Create a vector generator with given parameters.
+
+ Parameters
+ ----------
+ value : tuple[str, dict[str, Any]]
+ Class name and dictionary of input parameters
- :param configuration: the current configuration.
- :type configuration: a MDANSE.Framework.Configurable.Configurable object
- :param value: the Q vectors generator definition. It can be a 2-tuple, whose 1st element is the name of the Q vector generator \
- and 2nd element the parameters for this configurator or a string that matches a Q vectors user definition.
- :type value: 2-tuple or str
"""
self._original_input = value
trajConfig = self._configurable[self._dependencies["trajectory"]]
if isinstance(value, tuple):
try:
- generator, parameters = value
+ generator_name, parameters = value
except ValueError:
self.error_status = f"Invalid q vectors settings {value}"
return
generator = IQVectors.create(
- generator, trajConfig["instance"].configuration(0)
+ generator_name, trajConfig["instance"].configuration(0)
)
try:
generator.setup(parameters)
@@ -84,7 +83,7 @@ def configure(self, value):
if "q_vectors" not in generator.configuration:
self.error_status = "Wrong inputs for q-vector generation. At the moment there are no valid Q points."
return
- elif not generator.configuration["q_vectors"]:
+ if not generator.configuration["q_vectors"]:
self.error_status = "no Q vectors could be generated"
return
@@ -99,23 +98,33 @@ def configure(self, value):
self["shells"] = list(self["q_vectors"].keys())
self["n_shells"] = len(self["q_vectors"])
self["value"] = self["q_vectors"]
+ self["generator"] = generator
self.error_status = "OK"
def preview_output_axis(self):
+ """Output the values of |Q| from current parameters.
+
+ Returns
+ -------
+ list[float], str
+ Values of |Q| together with their physical unit
+
+ """
if not self.is_configured():
return None, None
if not self._valid:
return None, None
return self["shells"], "1/nm"
- def get_information(self):
- """
- Returns string information about this configurator.
+ def get_information(self) -> str:
+ """Return human-readable information about the generated vectors.
- :return: the information about this configurator.
- :rtype: str
- """
+ Returns
+ -------
+ str
+ Summary of generated vectors.
+ """
try:
info = [f"{self['n_shells']} Q shells generated\n"]
except KeyError:
diff --git a/MDANSE/Src/MDANSE/Framework/Jobs/CurrentCorrelationFunction.py b/MDANSE/Src/MDANSE/Framework/Jobs/CurrentCorrelationFunction.py
index 41aaa4bc2d..26d3e537e5 100644
--- a/MDANSE/Src/MDANSE/Framework/Jobs/CurrentCorrelationFunction.py
+++ b/MDANSE/Src/MDANSE/Framework/Jobs/CurrentCorrelationFunction.py
@@ -370,6 +370,10 @@ def finalize(self):
"""Normalize, Fourier transform and write the results out to
the MDA files.
"""
+ self.configuration["q_vectors"]["generator"].write_vectors_to_file(
+ self._outputData
+ )
+
nAtomsPerElement = self.configuration["atom_selection"].get_natoms()
for pair in self._elementsPairs:
pair_str = "".join(map(str, pair))
diff --git a/MDANSE/Src/MDANSE/Framework/Jobs/DynamicCoherentStructureFactor.py b/MDANSE/Src/MDANSE/Framework/Jobs/DynamicCoherentStructureFactor.py
index e04b575dcd..fedac5ff54 100644
--- a/MDANSE/Src/MDANSE/Framework/Jobs/DynamicCoherentStructureFactor.py
+++ b/MDANSE/Src/MDANSE/Framework/Jobs/DynamicCoherentStructureFactor.py
@@ -20,10 +20,14 @@
import numpy as np
from scipy.signal import correlate
+from MDANSE.MLogging import LOG
from MDANSE.Core.Error import Error
from MDANSE.Framework.Jobs.IJob import IJob
+from MDANSE.Framework.QVectors.LatticeQVectors import LatticeQVectors
from MDANSE.Mathematics.Arithmetic import assign_weights, get_weights, weighted_sum
from MDANSE.Mathematics.Signal import get_spectrum
+from MDANSE.Framework.QVectors.IQVectors import IQVectors
+from MDANSE.MolecularDynamics.UnitCell import UnitCell
class DynamicCoherentStructureFactorError(Error):
@@ -99,6 +103,12 @@ def initialize(self):
self.numberOfSteps = self.configuration["q_vectors"]["n_shells"]
nQShells = self.configuration["q_vectors"]["n_shells"]
+ if not isinstance(
+ self.configuration["q_vectors"]["generator"], LatticeQVectors
+ ):
+ LOG.warning(
+ f"This task should be used with a lattice-based Q vector generator. You have picked {self.configuration['q_vectors']['generator'].__class__}. The results are likely to be incorrect."
+ )
self._nFrames = self.configuration["frames"]["n_frames"]
@@ -182,6 +192,28 @@ def initialize(self):
main_result=True,
)
+ self._cell_std = 0.0
+ try:
+ all_cells = [
+ self.configuration["trajectory"]["instance"].unit_cell(frame)._unit_cell
+ for frame in self.configuration["frames"]["value"]
+ ]
+ except TypeError:
+ self._average_unit_cell = None
+ else:
+ self._average_unit_cell = UnitCell(
+ np.mean(
+ all_cells,
+ axis=0,
+ )
+ )
+ self._cell_std = UnitCell(
+ np.std(
+ all_cells,
+ axis=0,
+ )
+ )
+
def run_step(self, index):
"""
Runs a single step of the job.\n
@@ -212,9 +244,35 @@ def run_step(self, index):
dtype=np.complex64,
)
+ cell_present = True
+ cell_fixed = True
# loop over the trajectory time steps
for i, frame in enumerate(self.configuration["frames"]["value"]):
- qVectors = self.configuration["q_vectors"]["value"][shell]["q_vectors"]
+ unit_cell = traj.unit_cell(frame)
+ if unit_cell is None:
+ cell_present = False
+ elif not np.allclose(
+ unit_cell._unit_cell, self._average_unit_cell._unit_cell
+ ):
+ cell_fixed = False
+ if not cell_present:
+ qVectors = self.configuration["q_vectors"]["value"][shell][
+ "q_vectors"
+ ]
+ else:
+ try:
+ hkls = self.configuration["q_vectors"]["value"][shell]["hkls"]
+ except KeyError:
+ qVectors = self.configuration["q_vectors"]["value"][shell][
+ "q_vectors"
+ ]
+ else:
+ if hkls is None:
+ qVectors = self.configuration["q_vectors"]["value"][shell][
+ "q_vectors"
+ ]
+ else:
+ qVectors = IQVectors.hkl_to_qvectors(hkls, unit_cell)
coords = traj.configuration(frame)["coordinates"]
@@ -223,6 +281,14 @@ def run_step(self, index):
rho[element][i, :] = np.sum(
np.exp(1j * np.dot(selectedCoordinates, qVectors)), axis=0
)
+ if not cell_present:
+ LOG.warning(
+ "You are running the DCSF calculation on a trajectory without periodic boundary conditions."
+ )
+ if not cell_fixed:
+ LOG.warning(
+ f"The unit cell is VARIABLE with the standard deviation of {self._cell_std}. This analysis should not be used with NPT runs! PLEASE CHECK YOUR RESULTS CAREFULLY."
+ )
return index, rho
@@ -249,6 +315,10 @@ def finalize(self):
"""
Finalizes the calculations (e.g. averaging the total term, output files creations ...)
"""
+ self.configuration["q_vectors"]["generator"].write_vectors_to_file(
+ self._outputData
+ )
+
nAtomsPerElement = self.configuration["atom_selection"].get_natoms()
weights = self.configuration["weights"].get_weights()
weight_dict = get_weights(weights, nAtomsPerElement, 2, conc_exp=0.5)
diff --git a/MDANSE/Src/MDANSE/Framework/Jobs/DynamicIncoherentStructureFactor.py b/MDANSE/Src/MDANSE/Framework/Jobs/DynamicIncoherentStructureFactor.py
index 59fef14b6b..b1fb9ca25d 100644
--- a/MDANSE/Src/MDANSE/Framework/Jobs/DynamicIncoherentStructureFactor.py
+++ b/MDANSE/Src/MDANSE/Framework/Jobs/DynamicIncoherentStructureFactor.py
@@ -251,6 +251,9 @@ def finalize(self):
"""
Finalizes the calculations (e.g. averaging the total term, output files creations ...)
"""
+ self.configuration["q_vectors"]["generator"].write_vectors_to_file(
+ self._outputData
+ )
nAtomsPerElement = self.configuration["atom_selection"].get_natoms()
weights = self.configuration["weights"].get_weights()
diff --git a/MDANSE/Src/MDANSE/Framework/Jobs/ElasticIncoherentStructureFactor.py b/MDANSE/Src/MDANSE/Framework/Jobs/ElasticIncoherentStructureFactor.py
index 5f42bb6e17..3e3a3d602a 100644
--- a/MDANSE/Src/MDANSE/Framework/Jobs/ElasticIncoherentStructureFactor.py
+++ b/MDANSE/Src/MDANSE/Framework/Jobs/ElasticIncoherentStructureFactor.py
@@ -194,6 +194,10 @@ def finalize(self):
Finalizes the calculations (e.g. averaging the total term, output files creations ...)
"""
+ self.configuration["q_vectors"]["generator"].write_vectors_to_file(
+ self._outputData
+ )
+
nAtomsPerElement = self.configuration["atom_selection"].get_natoms()
for element, number in list(nAtomsPerElement.items()):
self._outputData[f"eisf_{element}"][:] /= number
diff --git a/MDANSE/Src/MDANSE/Framework/QVectors/ApproximateDispersionQVectors.py b/MDANSE/Src/MDANSE/Framework/QVectors/ApproximateDispersionQVectors.py
index 891203e473..de820e4bb1 100644
--- a/MDANSE/Src/MDANSE/Framework/QVectors/ApproximateDispersionQVectors.py
+++ b/MDANSE/Src/MDANSE/Framework/QVectors/ApproximateDispersionQVectors.py
@@ -24,7 +24,7 @@
class ApproximateDispersionQVectors(LatticeQVectors):
- """ """
+ """Generates Q vectors along a direction."""
settings = collections.OrderedDict()
settings["q_start"] = (
@@ -69,7 +69,7 @@ def _generate(self):
np.array(qStart)[:, np.newaxis] + np.outer(n, np.arange(0, nSteps)) * qStep
)
- hkls = np.rint(np.dot(self._directUnitCell, vects))
+ hkls = self.qvectors_to_hkl(vects, self._unit_cell)
dists = np.sqrt(np.sum(vects**2, axis=0))
dists = list(zip(range(len(dists)), dists))
@@ -94,5 +94,4 @@ def _generate(self):
if self._status is not None:
if self._status.is_stopped():
return
- else:
- self._status.update()
+ self._status.update()
diff --git a/MDANSE/Src/MDANSE/Framework/QVectors/CircularLatticeQVectors.py b/MDANSE/Src/MDANSE/Framework/QVectors/CircularLatticeQVectors.py
index 01175ad326..7e3e419143 100644
--- a/MDANSE/Src/MDANSE/Framework/QVectors/CircularLatticeQVectors.py
+++ b/MDANSE/Src/MDANSE/Framework/QVectors/CircularLatticeQVectors.py
@@ -19,13 +19,21 @@
import numpy as np
-from MDANSE.Mathematics.LinearAlgebra import Vector
-
from MDANSE.Framework.QVectors.LatticeQVectors import LatticeQVectors
+from MDANSE.Mathematics.LinearAlgebra import Vector
class CircularLatticeQVectors(LatticeQVectors):
- """ """
+ """Generates Q vectors on a plane.
+
+ Vectors are grouped into annuli (called 'shells')
+ based on their length.
+
+ Only vectors commensurate with the reciprocal
+ space lattice vectors will be generated.
+ |Q| values for which no valid vectors can
+ be found are omitted in the output.
+ """
settings = collections.OrderedDict()
settings["seed"] = ("IntegerConfigurator", {"mini": 0, "default": 0})
@@ -61,7 +69,7 @@ def _generate(self):
]
)
- qVects = np.dot(self._inverseUnitCell, hkls)
+ qVects = self.hkl_to_qvectors(hkls, self._unit_cell)
qMax = (
self._configuration["shells"]["last"]
@@ -109,15 +117,11 @@ def _generate(self):
self._configuration["q_vectors"][q]["q_vectors"] = vects[:, hits]
self._configuration["q_vectors"][q]["n_q_vectors"] = n
self._configuration["q_vectors"][q]["q"] = q
- self._configuration["q_vectors"][q]["hkls"] = np.rint(
- np.dot(
- self._directUnitCell,
- self._configuration["q_vectors"][q]["q_vectors"],
- )
+ self._configuration["q_vectors"][q]["hkls"] = self.qvectors_to_hkl(
+ vects[:, hits], self._unit_cell
)
if self._status is not None:
if self._status.is_stopped():
return
- else:
- self._status.update()
+ self._status.update()
diff --git a/MDANSE/Src/MDANSE/Framework/QVectors/CircularQVectors.py b/MDANSE/Src/MDANSE/Framework/QVectors/CircularQVectors.py
index b8a2b25a7e..67e727a25f 100644
--- a/MDANSE/Src/MDANSE/Framework/QVectors/CircularQVectors.py
+++ b/MDANSE/Src/MDANSE/Framework/QVectors/CircularQVectors.py
@@ -18,12 +18,16 @@
import numpy as np
-from MDANSE.Framework.QVectors.IQVectors import IQVectors, QVectorsError
+from MDANSE.Framework.QVectors.IQVectors import IQVectors
from MDANSE.Mathematics.Geometry import random_points_on_circle
class CircularQVectors(IQVectors):
- """ """
+ """Generates Q vectors on a plane.
+
+ Vectors are grouped into annuli (called 'shells')
+ based on their length.
+ """
settings = collections.OrderedDict()
settings["seed"] = ("IntegerConfigurator", {"mini": 0, "default": 0})
@@ -51,14 +55,11 @@ def _generate(self):
if self._configuration["seed"]["value"] != 0:
np.random.seed(self._configuration["seed"]["value"])
- try:
- axis = (
- self._configuration["axis_1"]["vector"]
- .cross(self._configuration["axis_2"]["vector"])
- .normal()
- )
- except ZeroDivisionError as e:
- raise QVectorsError(str(e))
+ axis = (
+ self._configuration["axis_1"]["vector"]
+ .cross(self._configuration["axis_2"]["vector"])
+ .normal()
+ )
width = self._configuration["width"]["value"]
@@ -84,5 +85,4 @@ def _generate(self):
if self._status is not None:
if self._status.is_stopped():
return
- else:
- self._status.update()
+ self._status.update()
diff --git a/MDANSE/Src/MDANSE/Framework/QVectors/DispersionLatticeQVectors.py b/MDANSE/Src/MDANSE/Framework/QVectors/DispersionLatticeQVectors.py
index ef4cc82071..48d126c391 100644
--- a/MDANSE/Src/MDANSE/Framework/QVectors/DispersionLatticeQVectors.py
+++ b/MDANSE/Src/MDANSE/Framework/QVectors/DispersionLatticeQVectors.py
@@ -22,7 +22,7 @@
class DispersionLatticeQVectors(LatticeQVectors):
- """ """
+ """Generates Q vectors along a direction."""
settings = collections.OrderedDict()
settings["start"] = (
@@ -48,7 +48,7 @@ def _generate(self):
)
# The k matrix (3,n_hkls)
- vects = np.dot(self._inverseUnitCell, hkls)
+ vects = self.hkl_to_qvectors(hkls, self._unit_cell)
dists = np.sqrt(np.sum(vects**2, axis=0))
@@ -69,5 +69,4 @@ def _generate(self):
if self._status is not None:
if self._status.is_stopped():
return
- else:
- self._status.update()
+ self._status.update()
diff --git a/MDANSE/Src/MDANSE/Framework/QVectors/GridQVectors.py b/MDANSE/Src/MDANSE/Framework/QVectors/GridQVectors.py
index 940cf2deec..c62c15e131 100644
--- a/MDANSE/Src/MDANSE/Framework/QVectors/GridQVectors.py
+++ b/MDANSE/Src/MDANSE/Framework/QVectors/GridQVectors.py
@@ -24,7 +24,14 @@
class GridQVectors(LatticeQVectors):
- """ """
+ """Generates vectors on a grid.
+
+ Vectors are generated from HKL values based on
+ the definition of the unit cell.
+
+ No symmetry considerations are used when generating
+ the vectors.
+ """
settings = collections.OrderedDict()
settings["hrange"] = (
@@ -59,7 +66,7 @@ def _generate(self):
hkls = hkls.reshape(3, nh * nk * nl)
# The k matrix (3,n_hkls)
- vects = np.dot(self._inverseUnitCell, hkls)
+ vects = self.hkl_to_qvectors(hkls, self._unit_cell)
dists = np.sqrt(np.sum(vects**2, axis=0))
@@ -75,7 +82,7 @@ def _generate(self):
dists.sort(key=operator.itemgetter(1))
qGroups = itertools.groupby(dists, key=operator.itemgetter(1))
qGroups = collections.OrderedDict(
- [(k, [item[0] for item in v]) for k, v in qGroups]
+ [(k, [item[0] for item in v]) for k, v in qGroups],
)
if self._status is not None:
@@ -93,5 +100,4 @@ def _generate(self):
if self._status is not None:
if self._status.is_stopped():
return
- else:
- self._status.update()
+ self._status.update()
diff --git a/MDANSE/Src/MDANSE/Framework/QVectors/IQVectors.py b/MDANSE/Src/MDANSE/Framework/QVectors/IQVectors.py
index b3fadad0ee..b77b230f27 100644
--- a/MDANSE/Src/MDANSE/Framework/QVectors/IQVectors.py
+++ b/MDANSE/Src/MDANSE/Framework/QVectors/IQVectors.py
@@ -14,18 +14,23 @@
# along with this program. If not, see .
#
import abc
+from typing import TYPE_CHECKING
+
+import numpy as np
from MDANSE.Core.Error import Error
-from MDANSE.Framework.Configurable import Configurable
from MDANSE.Core.SubclassFactory import SubclassFactory
+from MDANSE.Framework.Configurable import Configurable
from MDANSE.MLogging import LOG
-
-class QVectorsError(Error):
- pass
+if TYPE_CHECKING:
+ from MDANSE.Framework.OutputVariables.IOutputVariable import OutputData
+ from MDANSE.MolecularDynamics.UnitCell import UnitCell
class IQVectors(Configurable, metaclass=SubclassFactory):
+ """Parent class of all Q vector generators."""
+
is_lattice = False
def __init__(self, atom_configuration, status=None):
@@ -40,17 +45,110 @@ def _generate(self):
pass
def generate(self) -> bool:
+ """Generate vectors by calling the internal method _generate."""
if self._configured:
self._generate()
if self._status is not None:
self._status.finish()
return True
- else:
- LOG.error(
- "Cannot generate vectors: q vector generator is not configured correctly."
- )
- return False
+ LOG.error(
+ "Cannot generate vectors: q vector generator is not configured correctly.",
+ )
+ return False
- def setStatus(self, status):
- self._status = status
+ @classmethod
+ def qvectors_to_hkl(
+ self,
+ vector_array: np.array,
+ unit_cell: "UnitCell",
+ ) -> np.ndarray:
+ """Recalculate Q vectors to HKL Miller indices.
+
+ Using a unit cell definition, recalculates an array
+ of q vectors to an equivalent array of HKL Miller indices.
+
+ Parameters
+ ----------
+ vector_array : np.array
+ a (3,N) array of scattering vectors
+ unit_cell : UnitCell
+ an instance of UnitCell class describing the simulation box
+
+ Returns
+ -------
+ np.ndarray
+ A (3,N) array of HKL values (Miller indices)
+
+ """
+ return np.dot(unit_cell.direct, vector_array) / (2 * np.pi)
+
+ @classmethod
+ def hkl_to_qvectors(self, hkls: np.array, unit_cell: "UnitCell") -> np.ndarray:
+ """Convert an array of HKL values to scattering vectors.
+
+ Uses a unit cell object to get the lattice vectors for conversion.
+
+ Parameters
+ ----------
+ hkls : np.array
+ A (3,N) array of HKL values (Miller indices)
+ unit_cell : UnitCell
+ An instance of UnitCell class describing the simulation box shape
+
+ Returns
+ -------
+ np.ndarray
+ a (3, N) array of Q vectors (scattering vectors)
+
+ """
+ return 2 * np.pi * np.dot(unit_cell.inverse, hkls)
+
+ def write_vectors_to_file(self, output_data: "OutputData"):
+ """Write the vectors to output file as an array.
+
+ Writes a summary of the generated vectors to the output
+ file using an OutputData class instance.
+
+ Parameters
+ ----------
+ output_data : OutputData
+ An object managing the writeout to one or many output files
+
+ """
+ qvector_info = self._configuration["q_vectors"]
+ q_values = [float(x) for x in qvector_info]
+ output_data.add(
+ "vector_generator_q",
+ "LineOutputVariable",
+ q_values,
+ units="1/nm",
+ )
+ qvector_lengths = [qvector_info[q]["q_vectors"].shape[1] for q in q_values]
+ qarray_maxlength = np.max(qvector_lengths)
+ output_data.add(
+ "vector_generator_qvector_array",
+ "VolumeOutputVariable",
+ (len(q_values), 3, qarray_maxlength),
+ units="1/nm",
+ )
+ output_data["vector_generator_qvector_array"][:] = 0.0
+ for nq, q in enumerate(q_values):
+ output_data["vector_generator_qvector_array"][
+ nq, :, : qvector_lengths[nq]
+ ] = qvector_info[q]["q_vectors"]
+ try:
+ hkl_arrays = [qvector_info[q]["hkls"] for q in q_values]
+ except KeyError:
+ return
+ output_data.add(
+ "vector_generator_hkl_array",
+ "VolumeOutputVariable",
+ (len(q_values), 3, qarray_maxlength),
+ units="au",
+ )
+ output_data["vector_generator_hkl_array"][:] = 0.0
+ for nq, _ in enumerate(q_values):
+ output_data["vector_generator_hkl_array"][nq, :, : qvector_lengths[nq]] = (
+ hkl_arrays[nq]
+ )
diff --git a/MDANSE/Src/MDANSE/Framework/QVectors/LatticeQVectors.py b/MDANSE/Src/MDANSE/Framework/QVectors/LatticeQVectors.py
index 537f545156..985e643d84 100644
--- a/MDANSE/Src/MDANSE/Framework/QVectors/LatticeQVectors.py
+++ b/MDANSE/Src/MDANSE/Framework/QVectors/LatticeQVectors.py
@@ -14,25 +14,23 @@
# along with this program. If not, see .
#
-import numpy as np
-
-from MDANSE.Framework.QVectors.IQVectors import IQVectors, QVectorsError
+from MDANSE.Framework.QVectors.IQVectors import IQVectors
class LatticeQVectors(IQVectors):
+ """Parent class for vector generators which need unit cell information."""
+
is_lattice = True
def __init__(self, atom_configuration, status=None):
- super(LatticeQVectors, self).__init__(atom_configuration, status)
+ super().__init__(atom_configuration, status)
if atom_configuration is None:
- raise QVectorsError("No configuration set for the chemical system")
+ raise ValueError("No configuration set for the chemical system")
if not atom_configuration.is_periodic:
- raise QVectorsError(
+ raise ValueError(
"The universe must be periodic for building lattice-based Q vectors"
)
- self._inverseUnitCell = 2.0 * np.pi * atom_configuration.unit_cell.inverse
-
- self._directUnitCell = 2.0 * np.pi * atom_configuration.unit_cell.direct
+ self._unit_cell = atom_configuration.unit_cell
diff --git a/MDANSE/Src/MDANSE/Framework/QVectors/LinearLatticeQVectors.py b/MDANSE/Src/MDANSE/Framework/QVectors/LinearLatticeQVectors.py
index fefabc3aa3..26b2f43951 100644
--- a/MDANSE/Src/MDANSE/Framework/QVectors/LinearLatticeQVectors.py
+++ b/MDANSE/Src/MDANSE/Framework/QVectors/LinearLatticeQVectors.py
@@ -19,13 +19,25 @@
import numpy as np
-from MDANSE.Mathematics.LinearAlgebra import Vector
-
from MDANSE.Framework.QVectors.LatticeQVectors import LatticeQVectors
+from MDANSE.Mathematics.LinearAlgebra import Vector
class LinearLatticeQVectors(LatticeQVectors):
- """ """
+ """Generates vectors randomly on a straight line.
+
+ Only vectors commensurate with the reciprocal
+ space lattice vectors will be generated.
+ |Q| values for which no valid vectors can
+ be found are omitted in the output.
+
+ Vectors within one shell are generated within
+ a tolerance limit around a central |Q| value.
+ Most calculations will produce one data point
+ for |Q| by averaging the results over all
+ vectors in the group, which is still called
+ a shell.
+ """
settings = collections.OrderedDict()
settings["seed"] = ("IntegerConfigurator", {"mini": 0, "default": 0})
@@ -51,7 +63,9 @@ def _generate(self):
random.seed(self._configuration["seed"]["value"])
# The Q vector corresponding to the input hkl.
- qVect = np.dot(self._inverseUnitCell, self._configuration["axis"]["vector"])
+ qVect = self.hkl_to_qvectors(
+ self._configuration["axis"]["vector"], self._unit_cell
+ )
qMax = (
self._configuration["shells"]["last"]
@@ -95,15 +109,11 @@ def _generate(self):
self._configuration["q_vectors"][q]["q_vectors"] = vects[:, hits]
self._configuration["q_vectors"][q]["n_q_vectors"] = n
self._configuration["q_vectors"][q]["q"] = q
- self._configuration["q_vectors"][q]["hkls"] = np.rint(
- np.dot(
- self._directUnitCell,
- self._configuration["q_vectors"][q]["q_vectors"],
- )
+ self._configuration["q_vectors"][q]["hkls"] = self.qvectors_to_hkl(
+ vects[:, hits], self._unit_cell
)
if self._status is not None:
if self._status.is_stopped():
return
- else:
- self._status.update()
+ self._status.update()
diff --git a/MDANSE/Src/MDANSE/Framework/QVectors/LinearQVectors.py b/MDANSE/Src/MDANSE/Framework/QVectors/LinearQVectors.py
index 26c55932a1..26165ea4d3 100644
--- a/MDANSE/Src/MDANSE/Framework/QVectors/LinearQVectors.py
+++ b/MDANSE/Src/MDANSE/Framework/QVectors/LinearQVectors.py
@@ -22,7 +22,15 @@
class LinearQVectors(IQVectors):
- """ """
+ """Generates vectors randomly on a straight line.
+
+ Vectors within one shell are generated within
+ a tolerance limit around a central |Q| value.
+ Most calculations will produce one data point
+ for |Q| by averaging the results over all
+ vectors in the group, which is still called
+ a shell.
+ """
settings = collections.OrderedDict()
settings["seed"] = ("IntegerConfigurator", {"mini": 0, "default": 0})
@@ -73,5 +81,4 @@ def _generate(self):
if self._status is not None:
if self._status.is_stopped():
return
- else:
- self._status.update()
+ self._status.update()
diff --git a/MDANSE/Src/MDANSE/Framework/QVectors/MillerIndicesQVectors.py b/MDANSE/Src/MDANSE/Framework/QVectors/MillerIndicesQVectors.py
index 641b97bb24..d3c2c8e951 100644
--- a/MDANSE/Src/MDANSE/Framework/QVectors/MillerIndicesQVectors.py
+++ b/MDANSE/Src/MDANSE/Framework/QVectors/MillerIndicesQVectors.py
@@ -22,7 +22,13 @@
class MillerIndicesQVectors(LatticeQVectors):
- """ """
+ """Generates vectors on a grid.
+
+ Vectors are generated from HKL values based on
+ the definition of the unit cell.
+ They are then grouped into shells based on
+ their length.
+ """
settings = collections.OrderedDict()
settings["shells"] = (
@@ -70,7 +76,7 @@ def _generate(self):
hkls = hkls.reshape(3, int(round(hkls.size / 3)))
# The k matrix (3,n_hkls)
- vects = np.dot(self._inverseUnitCell, hkls)
+ vects = self.hkl_to_qvectors(hkls, self._unit_cell)
dists2 = np.sum(vects**2, axis=0)
@@ -96,15 +102,11 @@ def _generate(self):
self._configuration["q_vectors"][q]["q_vectors"] = vects[:, hits]
self._configuration["q_vectors"][q]["n_q_vectors"] = nHits
self._configuration["q_vectors"][q]["q"] = q
- self._configuration["q_vectors"][q]["hkls"] = np.rint(
- np.dot(
- self._directUnitCell,
- self._configuration["q_vectors"][q]["q_vectors"],
- )
+ self._configuration["q_vectors"][q]["hkls"] = self.qvectors_to_hkl(
+ vects[:, hits], self._unit_cell
)
if self._status is not None:
if self._status.is_stopped():
return
- else:
- self._status.update()
+ self._status.update()
diff --git a/MDANSE/Src/MDANSE/Framework/QVectors/SphericalLatticeQVectors.py b/MDANSE/Src/MDANSE/Framework/QVectors/SphericalLatticeQVectors.py
index 566170635a..0db2d05a5b 100644
--- a/MDANSE/Src/MDANSE/Framework/QVectors/SphericalLatticeQVectors.py
+++ b/MDANSE/Src/MDANSE/Framework/QVectors/SphericalLatticeQVectors.py
@@ -23,7 +23,19 @@
class SphericalLatticeQVectors(LatticeQVectors):
- """ """
+ """Generates vectors randomly on a sphere.
+
+ Only vectors commensurate with the reciprocal
+ space lattice vectors will be generated.
+ |Q| values for which no valid vectors can
+ be found are omitted in the output.
+
+ Vectors within one shell are generated within
+ a tolerance limit around a central |Q| value.
+ Most calculations will produce one data point
+ for |Q| by averaging the results over all
+ vectors in the shell.
+ """
settings = collections.OrderedDict()
settings["seed"] = ("IntegerConfigurator", {"mini": 0, "default": 0})
@@ -43,27 +55,25 @@ def _generate(self):
if self._configuration["seed"]["value"] != 0:
np.random.seed(self._configuration["seed"]["value"])
random.seed(self._configuration["seed"]["value"])
-
qMax = (
self._configuration["shells"]["last"]
+ 0.5 * self._configuration["width"]["value"]
)
- hklMax = (
- np.ceil([qMax / np.sqrt(np.sum(v**2)) for v in self._inverseUnitCell.T]) + 1
- )
+ hklMax = np.ceil(self.qvectors_to_hkl(qMax * np.eye(3), self._unit_cell)) + 1
- vects = np.mgrid[
- -hklMax[0] : hklMax[0] + 1,
- -hklMax[1] : hklMax[1] + 1,
- -hklMax[2] : hklMax[2] + 1,
+ hkl_vects = np.mgrid[
+ -hklMax[0, 0] : hklMax[0, 0] + 1,
+ -hklMax[1, 1] : hklMax[1, 1] + 1,
+ -hklMax[2, 2] : hklMax[2, 2] + 1,
]
- vects = vects.reshape(
- 3, int(2 * hklMax[0] + 1) * int(2 * hklMax[1] + 1) * int(2 * hklMax[2] + 1)
+ hkl_vects = hkl_vects.reshape(
+ 3,
+ np.prod(2 * np.diag(hklMax) + 1, dtype=int),
)
- vects = np.dot(self._inverseUnitCell, vects)
+ vects = self.hkl_to_qvectors(hkl_vects, self._unit_cell)
dists2 = np.sum(vects**2, axis=0)
@@ -96,15 +106,11 @@ def _generate(self):
self._configuration["q_vectors"][q]["q_vectors"] = vects[:, hits]
self._configuration["q_vectors"][q]["n_q_vectors"] = n
self._configuration["q_vectors"][q]["q"] = q
- self._configuration["q_vectors"][q]["hkls"] = np.rint(
- np.dot(
- self._directUnitCell,
- self._configuration["q_vectors"][q]["q_vectors"],
- )
+ self._configuration["q_vectors"][q]["hkls"] = self.qvectors_to_hkl(
+ vects[:, hits], self._unit_cell
)
if self._status is not None:
if self._status.is_stopped():
- return None
- else:
- self._status.update()
+ return
+ self._status.update()
diff --git a/MDANSE/Src/MDANSE/Framework/QVectors/SphericalQVectors.py b/MDANSE/Src/MDANSE/Framework/QVectors/SphericalQVectors.py
index d111e4072a..4478af08e4 100644
--- a/MDANSE/Src/MDANSE/Framework/QVectors/SphericalQVectors.py
+++ b/MDANSE/Src/MDANSE/Framework/QVectors/SphericalQVectors.py
@@ -18,12 +18,19 @@
import numpy as np
-from MDANSE.Mathematics.Geometry import random_points_on_sphere
from MDANSE.Framework.QVectors.IQVectors import IQVectors
+from MDANSE.Mathematics.Geometry import random_points_on_sphere
class SphericalQVectors(IQVectors):
- """ """
+ """Generates vectors randomly on a sphere.
+
+ Vectors within one shell are generated within
+ a tolerance limit around a central |Q| value.
+ Most calculations will produce one data point
+ for |Q| by averaging the results over all
+ vectors in the shell.
+ """
settings = collections.OrderedDict()
settings["seed"] = ("IntegerConfigurator", {"mini": 0, "default": 0})
@@ -68,5 +75,4 @@ def _generate(self):
if self._status is not None:
if self._status.is_stopped():
return
- else:
- self._status.update()
+ self._status.update()
diff --git a/MDANSE/Tests/UnitTests/Analysis/test_qvectors.py b/MDANSE/Tests/UnitTests/Analysis/test_qvectors.py
index 7c7f75117f..b51898281a 100644
--- a/MDANSE/Tests/UnitTests/Analysis/test_qvectors.py
+++ b/MDANSE/Tests/UnitTests/Analysis/test_qvectors.py
@@ -3,6 +3,8 @@
from os import path
import pytest
+import numpy as np
+
from MDANSE.Framework.InputData.HDFTrajectoryInputData import HDFTrajectoryInputData
from MDANSE.Framework.QVectors.IQVectors import IQVectors
from MDANSE.Framework.Jobs.IJob import IJob
@@ -22,6 +24,30 @@ def trajectory():
yield trajectory
+def test_qvector_to_hkl_conversion(trajectory):
+ for qvector_generator in IQVectors.indirect_subclasses():
+ instance = IQVectors.create(qvector_generator, trajectory.trajectory.configuration(0))
+ instance.setup({"shells": (5.0, 50.0, 10.0)})
+ unit_cell = trajectory.trajectory.unit_cell(0)
+ instance.generate()
+ try:
+ instance._configuration["shells"]
+ except KeyError:
+ print(f"{qvector_generator} has no shells")
+ continue
+ for q in instance._configuration["shells"]["value"][:2]:
+ try:
+ original_qvectors = instance._configuration["q_vectors"][q]["q_vectors"]
+ except KeyError:
+ continue
+ if len(original_qvectors) == 0:
+ continue
+ hkls = instance.qvectors_to_hkl(original_qvectors, unit_cell)
+ recalculated_qvectors = instance.hkl_to_qvectors(hkls, unit_cell)
+ assert np.allclose(original_qvectors, recalculated_qvectors)
+
+
+
def test_disf(trajectory):
temp_name = tempfile.mktemp()
parameters = {}
@@ -40,8 +66,6 @@ def test_disf(trajectory):
}
if len(qvector_defaults) < 1:
continue
- print(qvector_generator)
- print(qvector_defaults)
parameters["q_vectors"] = (qvector_generator, qvector_defaults)
disf = IJob.create("DynamicIncoherentStructureFactor")
disf.run(parameters, status=True)
diff --git a/MDANSE/Tests/UnitTests/Results/dcsf_short_traj.mda b/MDANSE/Tests/UnitTests/Results/dcsf_short_traj.mda
index 297a5a0405..18b73d0267 100644
Binary files a/MDANSE/Tests/UnitTests/Results/dcsf_short_traj.mda and b/MDANSE/Tests/UnitTests/Results/dcsf_short_traj.mda differ
diff --git a/MDANSE_GUI/Src/MDANSE_GUI/Resources/InstrumentDefinitions.toml b/MDANSE_GUI/Src/MDANSE_GUI/Resources/InstrumentDefinitions.toml
index 31755d5a5c..3253eef93e 100644
--- a/MDANSE_GUI/Src/MDANSE_GUI/Resources/InstrumentDefinitions.toml
+++ b/MDANSE_GUI/Src/MDANSE_GUI/Resources/InstrumentDefinitions.toml
@@ -5,7 +5,7 @@ _technique = "QENS"
_resolution_type = "Gaussian"
_resolution_fwhm = "0.1"
_resolution_unit = "meV"
-_qvector_type = "SphericalQVectors"
+_qvector_type = "SphericalLatticeQVectors"
_q_unit = "1/ang"
_q_min = "0.1"
_q_max = "4.0"
@@ -22,7 +22,7 @@ _technique = "QENS"
_resolution_type = "Gaussian"
_resolution_fwhm = "0.0254"
_resolution_unit = "meV"
-_qvector_type = "SphericalQVectors"
+_qvector_type = "SphericalLatticeQVectors"
_q_unit = "1/ang"
_q_min = "0.18"
_q_max = "1.8"
@@ -39,7 +39,7 @@ _technique = "QENS"
_resolution_type = "Gaussian"
_resolution_fwhm = "0.099"
_resolution_unit = "meV"
-_qvector_type = "SphericalQVectors"
+_qvector_type = "SphericalLatticeQVectors"
_q_unit = "1/ang"
_q_min = "0.37"
_q_max = "3.6"
@@ -56,7 +56,7 @@ _technique = "QENS"
_resolution_type = "Gaussian"
_resolution_fwhm = "0.0175"
_resolution_unit = "meV"
-_qvector_type = "SphericalQVectors"
+_qvector_type = "SphericalLatticeQVectors"
_q_unit = "1/ang"
_q_min = "0.18"
_q_max = "1.8"
@@ -73,7 +73,7 @@ _technique = "QENS"
_resolution_type = "Gaussian"
_resolution_fwhm = "0.0545"
_resolution_unit = "meV"
-_qvector_type = "SphericalQVectors"
+_qvector_type = "SphericalLatticeQVectors"
_q_unit = "1/ang"
_q_min = "0.37"
_q_max = "3.6"
@@ -90,7 +90,7 @@ _technique = "QENS"
_resolution_type = "ideal"
_resolution_fwhm = "0"
_resolution_unit = "meV"
-_qvector_type = "SphericalQVectors"
+_qvector_type = "SphericalLatticeQVectors"
_q_unit = "1/ang"
_q_min = "0.1"
_q_max = "10.0"
diff --git a/MDANSE_GUI/Src/MDANSE_GUI/Tabs/Models/PlottingContext.py b/MDANSE_GUI/Src/MDANSE_GUI/Tabs/Models/PlottingContext.py
index bd4d76f560..35dcb6c428 100644
--- a/MDANSE_GUI/Src/MDANSE_GUI/Tabs/Models/PlottingContext.py
+++ b/MDANSE_GUI/Src/MDANSE_GUI/Tabs/Models/PlottingContext.py
@@ -92,6 +92,11 @@ def __init__(self, name: str, source: "h5py.File"):
pass
self._axes = {}
self._axes_units = {}
+ if self._axes_tag == "index":
+ for dim_number, dim_length in enumerate(self._data.shape):
+ self._axes[f"index{dim_number}"] = np.arange(dim_length)
+ self._axes_units[f"index{dim_number}"] = "N/A"
+ return
for ax_number, axis_name in enumerate(self._axes_tag.split("|")):
aname = axis_name.strip()
if aname == "index":
diff --git a/MDANSE_GUI/Src/MDANSE_GUI/Tabs/Visualisers/Action.py b/MDANSE_GUI/Src/MDANSE_GUI/Tabs/Visualisers/Action.py
index a939f0eabe..1c507e26c8 100644
--- a/MDANSE_GUI/Src/MDANSE_GUI/Tabs/Visualisers/Action.py
+++ b/MDANSE_GUI/Src/MDANSE_GUI/Tabs/Visualisers/Action.py
@@ -335,7 +335,12 @@ def test_file_outputs(self):
def apply_instrument(self):
if self._current_instrument is not None:
- q_vector_tuple = self._current_instrument.create_q_vector_params()
+ initial_configuration = self._trajectory_configurator[
+ "instance"
+ ].configuration()
+ q_vector_tuple = self._current_instrument.create_q_vector_params(
+ initial_configuration
+ )
resolution_tuple = self._current_instrument.create_resolution_params()
for widget in self._widgets:
has_preview = callable(
diff --git a/MDANSE_GUI/Src/MDANSE_GUI/Tabs/Visualisers/InstrumentInfo.py b/MDANSE_GUI/Src/MDANSE_GUI/Tabs/Visualisers/InstrumentInfo.py
index 52738487d5..011cc7cdee 100644
--- a/MDANSE_GUI/Src/MDANSE_GUI/Tabs/Visualisers/InstrumentInfo.py
+++ b/MDANSE_GUI/Src/MDANSE_GUI/Tabs/Visualisers/InstrumentInfo.py
@@ -13,8 +13,9 @@
# You should have received a copy of the GNU General Public License
# along with this program. If not, see .
#
-from qtpy.QtCore import Signal, Slot, Qt
+
from qtpy.QtWidgets import QTextBrowser
+from qtpy.QtCore import Signal, Slot, Qt
from qtpy.QtGui import QStandardItem
from MDANSE.MLogging import LOG
@@ -82,7 +83,6 @@ def update_item(self):
def create_resolution_params(self):
if not self._configured:
return
- _general_resolution = ("gaussian", {"mu": 0.0, "sigma": 0.2})
calculator = ResolutionCalculator()
try:
calculator.update_model(self._resolution_type)
@@ -114,12 +114,12 @@ def sanitize_numbers(self):
results.append(new_entry)
return results
- def create_q_vector_params(self):
+ def create_q_vector_params(self, sample_configuration=None):
if not self._configured:
return
cov_type = self._qvector_type
try:
- qvec_generator = IQVectors.create(cov_type, None)
+ qvec_generator = IQVectors.create(cov_type, sample_configuration)
except ValueError:
return ("No qvectors", {})
except AttributeError:
@@ -162,7 +162,7 @@ def create_q_vector_params(self):
def filter_qvector_generator(self):
new_list = [str(x) for x in self.qvector_options]
if self._sample == "isotropic":
- new_list = [str(x) for x in new_list if "Lattice" not in x]
+ new_list = [str(x) for x in new_list if "Dispersion" not in x]
if self._technique == "QENS":
new_list = [str(x) for x in new_list if "Spherical" in x]
return new_list
@@ -170,7 +170,7 @@ def filter_qvector_generator(self):
def pick_qvector_generator(self) -> str:
qgenerator = "SphericalQVectors"
if self._sample == "isotropic" and self._technique == "QENS":
- qgenerator = "SphericalQVectors"
+ qgenerator = "SphericalLatticeQVectors"
return qgenerator