diff --git a/news/expose-aniso-pv.rst b/news/expose-aniso-pv.rst new file mode 100644 index 0000000..9e3f099 --- /dev/null +++ b/news/expose-aniso-pv.rst @@ -0,0 +1,25 @@ +**Added:** + +* Add explicitly constructable reflection profiles +* Expose objcryst's anisotropic Pseudo-Voigt profile +* Add PowderPatternDiffraction.SetProfile + +**Changed:** + +* + +**Deprecated:** + +* + +**Removed:** + +* + +**Fixed:** + +* + +**Security:** + +* diff --git a/src/extensions/powderpatterndiffraction_ext.cpp b/src/extensions/powderpatterndiffraction_ext.cpp index d29b822..def3afa 100644 --- a/src/extensions/powderpatterndiffraction_ext.cpp +++ b/src/extensions/powderpatterndiffraction_ext.cpp @@ -25,12 +25,24 @@ #include #include +#include #include namespace bp = boost::python; using namespace boost::python; using namespace ObjCryst; +namespace +{ + +void SetProfileCopy(PowderPatternDiffraction &pdiff, + const ReflectionProfile &profile) +{ + pdiff.SetProfile(profile.CreateCopy()); +} + +} // namespace + void wrap_powderpatterndiffraction() { @@ -60,6 +72,8 @@ void wrap_powderpatterndiffraction() (ReflectionProfile& (PowderPatternDiffraction::*)()) &PowderPatternDiffraction::GetProfile, return_internal_reference<>()) + .def("SetProfile", &SetProfileCopy, bp::arg("profile"), + "Install an independent copy of a reflection profile.") .def("SetExtractionMode", &PowderPatternDiffraction::SetExtractionMode, (bp::arg("extract")=true, bp::arg("init")=false)) diff --git a/src/extensions/reflectionprofile_ext.cpp b/src/extensions/reflectionprofile_ext.cpp index c734fbf..e72bc26 100644 --- a/src/extensions/reflectionprofile_ext.cpp +++ b/src/extensions/reflectionprofile_ext.cpp @@ -17,7 +17,10 @@ * *****************************************************************************/ +#define BOOST_PYTHON_MAX_ARITY 20 + #include +#include #include #include #undef B0 @@ -113,6 +116,8 @@ void wrap_reflectionprofile() (bp::arg("relativeIntensity"), bp::arg("xcenter"), bp::arg("h"), bp::arg("k"), bp::arg("l")), "Return the full profile width at a given relative intensity for reflection (h, k, l) around `xcenter`.") + .def("IsAnisotropic", &ReflectionProfile::IsAnisotropic, + "Return whether the profile depends on the reflection indices.") .def("XMLOutput", pure_virtual((void (ReflectionProfile::*)(ostream &, int) const) & ReflectionProfile::XMLOutput), (bp::arg("os"), bp::arg("indent")), @@ -121,4 +126,44 @@ void wrap_reflectionprofile() pure_virtual((void (ReflectionProfile::*)(istream &, const XMLCrystTag &))&ReflectionProfile::XMLInput), (bp::arg("is"), bp::arg("tag")), "Load ReflectionProfile parameters from an XML stream and tag."); + + class_ >( + "ReflectionProfilePseudoVoigt", init<>()) + .def(init(bp::arg("old"))) + .def("CreateCopy", &ReflectionProfilePseudoVoigt::CreateCopy, + return_value_policy(), + "Return an independent copy of this profile.") + .def("SetProfilePar", &ReflectionProfilePseudoVoigt::SetProfilePar, + (bp::arg("fwhmCagliotiW"), bp::arg("fwhmCagliotiU") = 0, + bp::arg("fwhmCagliotiV") = 0, bp::arg("eta0") = 0.5, + bp::arg("eta1") = 0), + "Set the isotropic pseudo-Voigt profile parameters."); + + class_ >( + "ReflectionProfilePseudoVoigtAnisotropic", init<>()) + .def(init( + bp::arg("old"))) + .def("CreateCopy", + &ReflectionProfilePseudoVoigtAnisotropic::CreateCopy, + return_value_policy(), + "Return an independent copy of this profile.") + .def( + "SetProfilePar", + &ReflectionProfilePseudoVoigtAnisotropic::SetProfilePar, + (bp::arg("fwhmCagliotiW"), bp::arg("fwhmCagliotiU") = 0, + bp::arg("fwhmCagliotiV") = 0, bp::arg("fwhmGaussP") = 0, + bp::arg("fwhmLorentzX") = 0, + bp::arg("fwhmLorentzY") = 0, + bp::arg("fwhmLorentzGammaHH") = 0, + bp::arg("fwhmLorentzGammaKK") = 0, + bp::arg("fwhmLorentzGammaLL") = 0, + bp::arg("fwhmLorentzGammaHK") = 0, + bp::arg("fwhmLorentzGammaHL") = 0, + bp::arg("fwhmLorentzGammaKL") = 0, + bp::arg("pseudoVoigtEta0") = 0, + bp::arg("pseudoVoigtEta1") = 0, + bp::arg("asymA0") = 1, bp::arg("asymA1") = 0, + bp::arg("asymA2") = 0), + "Set the anisotropic pseudo-Voigt profile parameters."); } diff --git a/src/pyobjcryst/reflectionprofile.py b/src/pyobjcryst/reflectionprofile.py index c749b57..c88f4b6 100644 --- a/src/pyobjcryst/reflectionprofile.py +++ b/src/pyobjcryst/reflectionprofile.py @@ -20,6 +20,22 @@ ``XMLOutput`` / ``XMLInput`` and ``CreateCopy``. ``GetProfile`` accepts a Python sequence or numpy array for ``x``. +Concrete isotropic and anisotropic pseudo-Voigt profiles can be installed +with ``PowderPatternDiffraction.SetProfile``. The diffraction component +stores an independent copy, so one configured profile can be reused as a +template for multiple phases:: + + from pyobjcryst.reflectionprofile import ( + ReflectionProfilePseudoVoigtAnisotropic, + ) + + profile = ReflectionProfilePseudoVoigtAnisotropic() + profile.GetPar("W").SetValue(1e-6) + profile.GetPar("G_HH").SetValue(2e-6) + + for pdiff in powder_pattern.get_crystalline_components(): + pdiff.SetProfile(profile) + Example ------- @@ -48,6 +64,16 @@ y_broader = profile.GetProfile(window, xcenter, 1, 0, 0) """ -__all__ = ["ReflectionProfile", "ReflectionProfileType"] +__all__ = [ + "ReflectionProfile", + "ReflectionProfilePseudoVoigt", + "ReflectionProfilePseudoVoigtAnisotropic", + "ReflectionProfileType", +] -from pyobjcryst._pyobjcryst import ReflectionProfile, ReflectionProfileType +from pyobjcryst._pyobjcryst import ( + ReflectionProfile, + ReflectionProfilePseudoVoigt, + ReflectionProfilePseudoVoigtAnisotropic, + ReflectionProfileType, +) diff --git a/tests/test_reflectionprofile.py b/tests/test_reflectionprofile.py index a826e51..c304e44 100644 --- a/tests/test_reflectionprofile.py +++ b/tests/test_reflectionprofile.py @@ -7,6 +7,11 @@ from pyobjcryst.powderpattern import PowderPattern from pyobjcryst.refinableobj import RefinableObj +from pyobjcryst.reflectionprofile import ( + ReflectionProfile, + ReflectionProfilePseudoVoigt, + ReflectionProfilePseudoVoigtAnisotropic, +) class TestReflectionProfile(unittest.TestCase): @@ -19,17 +24,120 @@ def prepare_fixture(self, loadcifdata): def setUp(self): """Set up a ReflectionProfile instance for testing.""" x = np.linspace(0, 40, 1000) - c = self.loadcifdata("paracetamol.cif") + self.crystal = self.loadcifdata("paracetamol.cif") self.pp = PowderPattern() self.pp.SetWavelength(0.7) self.pp.SetPowderPatternX(np.deg2rad(x)) self.pp.SetPowderPatternObs(np.ones_like(x)) - self.ppd = self.pp.AddPowderPatternDiffraction(c) + self.ppd = self.pp.AddPowderPatternDiffraction(self.crystal) self.profile = self.ppd.GetProfile() + def test_concrete_pseudo_voigt_profiles(self): + """Concrete isotropic and anisotropic profiles are + constructible.""" + isotropic = ReflectionProfilePseudoVoigt() + anisotropic = ReflectionProfilePseudoVoigtAnisotropic() + + self.assertIsInstance(isotropic, ReflectionProfile) + self.assertIsInstance(anisotropic, ReflectionProfile) + self.assertFalse(isotropic.IsAnisotropic()) + self.assertTrue(anisotropic.IsAnisotropic()) + + expected = { + "U", + "V", + "W", + "P", + "X", + "Y", + "G_HH", + "G_KK", + "G_LL", + "G_HK", + "G_HL", + "G_KL", + "Eta0", + "Eta1", + "Asym0", + "Asym1", + "Asym2", + } + self.assertEqual( + expected, + { + anisotropic.GetPar(i).GetName() + for i in range(anisotropic.GetNbPar()) + }, + ) + + def test_anisotropic_set_profile_par(self): + """SetProfilePar assigns widths and retains symmetric default + asymmetry.""" + profile = ReflectionProfilePseudoVoigtAnisotropic() + profile.SetProfilePar(1e-6, fwhmLorentzX=2e-6) + profile.GetPar("P").SetValue(3e-6) + + self.assertAlmostEqual(profile.GetPar("W").GetValue(), 1e-6) + self.assertAlmostEqual(profile.GetPar("X").GetValue(), 2e-6) + self.assertAlmostEqual(profile.GetPar("P").GetValue(), 3e-6) + self.assertAlmostEqual(profile.GetPar("Asym0").GetValue(), 1.0) + + def test_anisotropic_profile_depends_on_hkl(self): + """Anisotropic Lorentz coefficients produce direction + dependence.""" + profile = ReflectionProfilePseudoVoigtAnisotropic() + profile.SetProfilePar( + 1e-6, + fwhmLorentzGammaHH=1e-3, + fwhmLorentzGammaKK=2e-3, + pseudoVoigtEta0=1, + ) + center = 0.5 + x = np.linspace(center - 0.1, center + 0.1, 401) + + h00 = profile.GetProfile(x, center, 1, 0, 0) + zero_k0 = profile.GetProfile(x, center, 0, 1, 0) + + self.assertFalse(np.allclose(h00, zero_k0)) + + def test_set_profile_copies_reusable_template(self): + """One profile template can safely initialize multiple + phases.""" + second = self.pp.AddPowderPatternDiffraction(self.crystal) + template = ReflectionProfilePseudoVoigtAnisotropic() + template.GetPar("W").SetValue(1e-6) + template.GetPar("G_HH").SetValue(2e-6) + + for pdiff in (self.ppd, second): + pdiff.SetProfile(template) + + first_profile = self.ppd.GetProfile() + second_profile = second.GetProfile() + self.assertIsInstance( + first_profile, ReflectionProfilePseudoVoigtAnisotropic + ) + self.assertIsInstance( + second_profile, ReflectionProfilePseudoVoigtAnisotropic + ) + self.assertAlmostEqual(first_profile.GetPar("W").GetValue(), 1e-6) + self.assertAlmostEqual(second_profile.GetPar("G_HH").GetValue(), 2e-6) + + first_profile.GetPar("W").SetValue(3e-6) + self.assertAlmostEqual(template.GetPar("W").GetValue(), 1e-6) + self.assertAlmostEqual(second_profile.GetPar("W").GetValue(), 1e-6) + + def test_set_profile_accepts_temporary(self): + """A temporary concrete profile is safely copied into the + phase.""" + self.ppd.SetProfile(ReflectionProfilePseudoVoigt()) + profile = self.ppd.GetProfile() + + self.assertIsInstance(profile, ReflectionProfilePseudoVoigt) + self.assertFalse(profile.IsAnisotropic()) + def test_get_computed_profile(self): """Sample a profile slice and verify broadening lowers the peak height."""