From 8fe31e25bbcc426f280c45db29dca8aac020fc66 Mon Sep 17 00:00:00 2001 From: Caden Myers Date: Tue, 27 Oct 2026 15:10:34 -0400 Subject: [PATCH] feat: expose March-Dollase preferred-orientation correction Bind TextureMarchDollase in powderpatterndiffraction_ext.cpp and add PowderPatternDiffraction.GetTextureMarchDollase, which finds the component's texture object in its sub-object registry because ObjCryst++ has no public accessor. SetPhasePar sets values through the phase's RefinablePars, since ObjCryst++'s SetPhasePar does not click the texture clock and leaves the cached pattern stale. DeletePhase is not bound because libObjCryst does not define it. Co-Authored-By: Claude Opus 5.5 --- news/march-dollase-texture.rst | 23 +++ .../powderpatterndiffraction_ext.cpp | 57 +++++++ src/pyobjcryst/powderpattern.py | 2 + tests/test_powderpattern.py | 146 +++++++++++++++++- 4 files changed, 227 insertions(+), 1 deletion(-) create mode 100644 news/march-dollase-texture.rst diff --git a/news/march-dollase-texture.rst b/news/march-dollase-texture.rst new file mode 100644 index 0000000..d62cecf --- /dev/null +++ b/news/march-dollase-texture.rst @@ -0,0 +1,23 @@ +**Added:** + +* Expose the March-Dollase preferred-orientation correction TextureMarchDollase and PowderPatternDiffraction.GetTextureMarchDollase. + +**Changed:** + +* + +**Deprecated:** + +* + +**Removed:** + +* + +**Fixed:** + +* + +**Security:** + +* diff --git a/src/extensions/powderpatterndiffraction_ext.cpp b/src/extensions/powderpatterndiffraction_ext.cpp index 839c925..d504e99 100644 --- a/src/extensions/powderpatterndiffraction_ext.cpp +++ b/src/extensions/powderpatterndiffraction_ext.cpp @@ -26,8 +26,11 @@ #include #include #include +#include #include +#include + namespace bp = boost::python; using namespace boost::python; using namespace ObjCryst; @@ -41,6 +44,40 @@ void SetProfileCopy(PowderPatternDiffraction &pdiff, pdiff.SetProfile(profile.CreateCopy()); } +// PowderPatternDiffraction has no public accessor for its texture +// correction, but registers it as a sub-object. +TextureMarchDollase &GetTextureMarchDollase(PowderPatternDiffraction &pdiff) +{ + ObjRegistry &subobjs = pdiff.GetSubObjRegistry(); + for (long i = 0; i < subobjs.GetNb(); ++i) + { + TextureMarchDollase *texture = + dynamic_cast(&subobjs.GetObj(i)); + if (texture != 0) return *texture; + } + throw ObjCrystException("PowderPatternDiffraction has no " + "TextureMarchDollase sub-object. Please report " + "this as a pyobjcryst bug."); +} + +// TextureMarchDollase::SetPhasePar does not click the texture clock, so the +// cached powder pattern is not recalculated. Set the values through the +// phase's RefinablePars instead, which AddPhase registers in this order. +void SetPhasePar(TextureMarchDollase &texture, const unsigned int i, + const REAL fraction, const REAL coeffMarch, + const REAL h, const REAL k, const REAL l) +{ + if (i >= texture.GetNbPhase()) + { + throw ObjCrystException("TextureMarchDollase has no phase " + + std::to_string(i) + ". Please use an index " + "smaller than GetNbPhase()."); + } + const REAL values[5] = {fraction, coeffMarch, h, k, l}; + for (long j = 0; j < 5; ++j) + texture.GetPar(5L * i + j).SetValue(values[j]); +} + } // namespace @@ -91,6 +128,26 @@ void wrap_powderpatterndiffraction() .def("GetFhklObsSq", &PowderPatternDiffraction::GetFhklObsSq, return_value_policy()) .def("X2XCorrPhase", &PowderPatternDiffraction::X2XCorrPhase) + .def("GetTextureMarchDollase", &GetTextureMarchDollase, + return_internal_reference<>(), + "Return the March-Dollase preferred-orientation correction.") + ; + + class_, boost::noncopyable>( + "TextureMarchDollase", no_init) + .def("AddPhase", &TextureMarchDollase::AddPhase, + (bp::arg("fraction"), bp::arg("coeffMarch"), bp::arg("h"), + bp::arg("k"), bp::arg("l"))) + .def("SetPhasePar", &SetPhasePar, + (bp::arg("i"), bp::arg("fraction"), bp::arg("coeffMarch"), + bp::arg("h"), bp::arg("k"), bp::arg("l"))) + .def("GetNbPhase", &TextureMarchDollase::GetNbPhase) + .def("GetFraction", &TextureMarchDollase::GetFraction, bp::arg("i")) + .def("GetMarchCoeff", &TextureMarchDollase::GetMarchCoeff, + bp::arg("i")) + .def("GetPhaseH", &TextureMarchDollase::GetPhaseH, bp::arg("i")) + .def("GetPhaseK", &TextureMarchDollase::GetPhaseK, bp::arg("i")) + .def("GetPhaseL", &TextureMarchDollase::GetPhaseL, bp::arg("i")) ; } diff --git a/src/pyobjcryst/powderpattern.py b/src/pyobjcryst/powderpattern.py index a044024..0b74b26 100644 --- a/src/pyobjcryst/powderpattern.py +++ b/src/pyobjcryst/powderpattern.py @@ -30,6 +30,7 @@ "PowderPatternDiffraction", "ReflectionProfileType", "SpaceGroupExplorer", + "TextureMarchDollase", ] from pyobjcryst._pyobjcryst import ( @@ -45,6 +46,7 @@ PowderPatternDiffraction, ReflectionProfileType, SpaceGroupExplorer, + TextureMarchDollase, ) from pyobjcryst.general import ObjCrystException from pyobjcryst.refinableobj import refpartype_scattdata_background diff --git a/tests/test_powderpattern.py b/tests/test_powderpattern.py index a5420a1..02ba692 100644 --- a/tests/test_powderpattern.py +++ b/tests/test_powderpattern.py @@ -14,6 +14,7 @@ ############################################################################## """Unit tests for pyobjcryst.powderpattern (with indexing &""" +import re import unittest from unittest.mock import MagicMock, patch @@ -23,7 +24,11 @@ from pyobjcryst import ObjCrystException from pyobjcryst.indexing import CrystalCentering, CrystalSystem, quick_index -from pyobjcryst.powderpattern import PowderPattern, SpaceGroupExplorer +from pyobjcryst.powderpattern import ( + PowderPattern, + SpaceGroupExplorer, + TextureMarchDollase, +) from pyobjcryst.radiation import RadiationType, WavelengthType from pyobjcryst.reflectionprofile import ReflectionProfileType @@ -408,5 +413,144 @@ def test_X2XCorrPhase_flat_detector_displacement(self): # ---------------------------------------------------------------------------- + +def build_textured_pattern(crystal, phases): + """Return an X-ray powder pattern of crystal with March-Dollase + phases.""" + pattern = PowderPattern() + pattern.SetWavelength(1.54) + x = np.deg2rad(np.linspace(5, 40, 2001)) + pattern.SetPowderPatternX(x) + pattern.SetPowderPatternObs(np.ones_like(x)) + diffraction = pattern.AddPowderPatternDiffraction(crystal) + texture = diffraction.GetTextureMarchDollase() + for phase in phases: + texture.AddPhase(*phase) + return pattern + + +@pytest.mark.parametrize( + "phases, expected_par_names", + [ + # Test that AddPhase stores a texture phase and its refinable + # parameters + # C1: One phase with a fraction, March coefficient and hkl is + # added. + # Expected: The getters return the phase's values, and its five + # parameters are named with the suffix _0. + ( + [(1.0, 0.5, 0, 0, 1)], + ["Fraction_0", "MarchCoeff_0", "H_0", "K_0", "L_0"], + ), + # C2: Two phases are added. + # Expected: The getters return each phase's values, and the second + # phase's parameters are named with the suffix _1. + ( + [(0.7, 0.5, 0, 0, 1), (0.3, 2.0, 1, 1, 0)], + [ + "Fraction_0", + "MarchCoeff_0", + "H_0", + "K_0", + "L_0", + "Fraction_1", + "MarchCoeff_1", + "H_1", + "K_1", + "L_1", + ], + ), + ], +) +def test_add_phase(loadcifdata, phases, expected_par_names): + crystal = loadcifdata("paracetamol.cif") + pattern = build_textured_pattern(crystal, phases) + texture = pattern.GetPowderPatternComponent(0).GetTextureMarchDollase() + actual_type = type(texture) + expected_type = TextureMarchDollase + actual_phases = [ + ( + texture.GetFraction(i), + texture.GetMarchCoeff(i), + texture.GetPhaseH(i), + texture.GetPhaseK(i), + texture.GetPhaseL(i), + ) + for i in range(texture.GetNbPhase()) + ] + expected_phases = phases + actual_par_names = [ + texture.GetPar(i).GetName() for i in range(texture.GetNbPar()) + ] + assert actual_type == expected_type + assert actual_phases == pytest.approx(expected_phases) + assert actual_par_names == expected_par_names + + +@pytest.mark.parametrize( + "initial_phases, method_name, method_args, expected_phases", + [ + # Test that changing the texture of an existing pattern + # recalculates it + # C1: A phase is added to an untextured pattern. + # Expected: The pattern equals one built with that phase. + ([], "AddPhase", (1.0, 0.5, 0, 0, 1), [(1.0, 0.5, 0, 0, 1)]), + # C2: SetPhasePar changes the March coefficient of a phase. + # Expected: The pattern equals one built with the new March + # coefficient. + ( + [(1.0, 0.5, 0, 0, 1)], + "SetPhasePar", + (0, 1.0, 2.0, 0, 0, 1), + [(1.0, 2.0, 0, 0, 1)], + ), + ], +) +def test_get_powder_pattern_calc_after_texture_change( + loadcifdata, initial_phases, method_name, method_args, expected_phases +): + crystal = loadcifdata("paracetamol.cif") + pattern = build_textured_pattern(crystal, initial_phases) + pattern.GetPowderPatternCalc() + texture = pattern.GetPowderPatternComponent(0).GetTextureMarchDollase() + getattr(texture, method_name)(*method_args) + actual_calc = pattern.GetPowderPatternCalc() + expected_calc = build_textured_pattern( + crystal, expected_phases + ).GetPowderPatternCalc() + assert np.allclose(actual_calc, expected_calc) + + +def test_set_phase_par_out_of_range(loadcifdata): + # Test SetPhasePar with an index of a phase that does not exist + # Expected: An ObjCrystException names the missing phase and says + # which indices are valid. + crystal = loadcifdata("paracetamol.cif") + pattern = build_textured_pattern(crystal, [(1.0, 0.5, 0, 0, 1)]) + texture = pattern.GetPowderPatternComponent(0).GetTextureMarchDollase() + expected_msg = ( + "TextureMarchDollase has no phase 1. Please use an index smaller " + "than GetNbPhase()." + ) + with pytest.raises(ObjCrystException, match=re.escape(expected_msg)): + texture.SetPhasePar(1, 1.0, 2.0, 0, 0, 1) + + +def test_get_powder_pattern_calc_after_march_coeff_par_change(loadcifdata): + # Test refining the March coefficient through its RefinablePar + # Expected: The pattern equals one built with the new March + # coefficient. + crystal = loadcifdata("paracetamol.cif") + pattern = build_textured_pattern(crystal, [(1.0, 0.5, 0, 0, 1)]) + pattern.GetPowderPatternCalc() + texture = pattern.GetPowderPatternComponent(0).GetTextureMarchDollase() + texture.GetPar("MarchCoeff_0").SetValue(2.0) + actual_calc = pattern.GetPowderPatternCalc() + expected_calc = build_textured_pattern( + crystal, [(1.0, 2.0, 0, 0, 1)] + ).GetPowderPatternCalc() + assert np.allclose(actual_calc, expected_calc) + + if __name__ == "__main__": unittest.main()