add "wavelength" to Diffractometer class

[prjemian]

As suggested, add wavelength to hkl.diffract.Diffractometer, comparable to the energy attribute. Consider deprecating the energy attribute since it is implemented specific to X-ray probes.

Since wavelength is the scientific term used in the computations, consider deprecating the energy attribute since energy is implemented here as specific to X-ray probes. If energy is to be deprecated, the final removal of energy from hklpy would be a breaking change (upgrade the major version number). That should be considered with other changes (such as the consideration of units handling for wavelength and lattice parameters).

[prjemian]

These are the places where energy is mentioned in the hkl c++ source code, master branch:

(base) prjemian@poof ~/.../projects/hkl-cpp $ git grep -i energy
Documentation/hkl.org.in:     =HklEngine= depends of the axes, the energy, or the sample. the
Documentation/hkl.org.in:  if (dependencies & HKL_ENGINE_DEPENDENCIES_ENERGY) {
hkl.h:  HKL_ENGINE_DEPENDENCIES_ENERGY = 1u << 1,
hkl/hkl-engine-petra3-p08-lisa.c:For a given energy, E, one can compute mchi, sphi, dtth, dh and drot
hkl/hkl-engine-petra3-p08-lisa.c:lambda = hc_over_e/Energy
hkl/hkl-engine-petra3-p08-lisa.c:the energy is changed, in order to correct the shifts in the position
hkl/hkl-pseudoaxis-common-hkl.c:                                HKL_ENGINE_DEPENDENCIES_AXES | HKL_ENGINE_DEPENDENCIES_ENERGY | HKL_ENGINE_DEPENDENCIES_SAMPLE),
hkl/hkl-pseudoaxis-common-psi.c:                                HKL_ENGINE_DEPENDENCIES_AXES | HKL_ENGINE_DEPENDENCIES_ENERGY | HKL_ENGINE_DEPENDENCIES_SAMPLE),
hkl/hkl-pseudoaxis-common-q.c:                          HKL_ENGINE_DEPENDENCIES_AXES | HKL_ENGINE_DEPENDENCIES_ENERGY),
hkl/hkl-pseudoaxis-common-q.c:                          HKL_ENGINE_DEPENDENCIES_AXES | HKL_ENGINE_DEPENDENCIES_ENERGY),
hkl/hkl-pseudoaxis-common-q.c:                          HKL_ENGINE_DEPENDENCIES_AXES | HKL_ENGINE_DEPENDENCIES_ENERGY),
tests/bindings/crystal.ini:AutoEnergyUpdate 1
tests/bindings/polarisation.py:Source = namedtuple('Source', ['wavelength', 'energy'])
tests/bindings/polarisation.py:                energy = 12.39842 / wavelength
tests/bindings/polarisation.py:        source = Source(wavelength, energy)
tests/bindings/polarisation.py:        print('Energy (keV):', config.source)
tests/hkl-pseudoaxis-t.c:               res &= DIAG((depends & (HKL_ENGINE_DEPENDENCIES_AXES | HKL_ENGINE_DEPENDENCIES_ENERGY | HKL_ENGINE_DEPENDENCIES_SAMPLE)) != 0);
tests/hkl-pseudoaxis-t.c:               res &= DIAG((depends & HKL_ENGINE_DEPENDENCIES_ENERGY) == 0);
tests/hkl-pseudoaxis-t.c:               res &= DIAG((depends &(HKL_ENGINE_DEPENDENCIES_AXES | HKL_ENGINE_DEPENDENCIES_ENERGY)) != 0);
tests/hkl-pseudoaxis-t.c:               res &= DIAG((depends &(HKL_ENGINE_DEPENDENCIES_AXES | HKL_ENGINE_DEPENDENCIES_ENERGY)) != 0);
tests/hkl-pseudoaxis-t.c:               res &= DIAG((depends & (HKL_ENGINE_DEPENDENCIES_AXES | HKL_ENGINE_DEPENDENCIES_ENERGY)) != 0);
tests/hkl-pseudoaxis-t.c:               res &= DIAG((depends & (HKL_ENGINE_DEPENDENCIES_AXES | HKL_ENGINE_DEPENDENCIES_ENERGY | HKL_ENGINE_DEPENDENCIES_SAMPLE)) != 0);

[prjemian]

Wavelength would be implemented similar to energy.

The energy attribute of class hkl.diffract.Diffractometer() is for convenience when integrating with the control system for a specific instrument. Locally, it is expected to subclass from Diffractometer and re-define energy, such as:

class LocalDiffractometer(hkl.diffract.Diffractometer):
    # ...
    energy = Component(EpicsSignal, "EPICS:diffractometer:energy")
    # ...

Then re-define this method: https://github.com/bluesky/hklpy/blob/9f24b5d025ee7a185a6c152554f0d4976e3da76c/hkl/diffract.py#L140-L146

and handle any conversion of the source units (EPICS:diffractometer:energy.EGU) to the keV as required by hkl.calc.energy.

[prjemian]

@gfabbris, @ambarb - Discussion?

[prjemian]

Of course, this means that PV EPICS:diffractometer:energy.EGU must have a units string that can be converted to keV.

[ambarb]

I need to understand a little better how we want to link the beamline energy to tardis.calc.energy (I think you used tardis.energy.put()?). Right now, we simply setting equal to a float. We also do not want to permanently link the beamline energy to the hklpy_object.energy

We also need to think about energy calibration. Maybe this requires a second PV. At CSX, it was decided that we would not adjust lattice constants, but we can offset the energy by adding the energy offset (usually up to +/-3 eV). tardis.calc.energy = N where N contains the energy offset, pgm.energy.user_setpoint, and our unit conversion factors.

[ambarb]

Also there seems to be 2 approaches. The concerns for approach 1 are above. For approach 2, it seems a little easier to manage and doesn't require any one to make new PVs.

1. Use EGU of beamline energy
2. Set unit configuration attributes that are hand configured by the beamline:
   • beamline energy
   • unit cell length and angle units

With both approaches we could still think about making a "soft" link between beamline energy and the hklpy object. This may help in getting a standard "fixedQenergy" object that can be "scanned" or a default plan that comes with hklpy. WE would have to import the RunEngine, but this isn't a bad idea either because we need to add RE integration tests instead of manually doing them at the beamline.

[prjemian]

I agree that 2 is the most popular convention now. An example could demonstrate 1 and include an energy offset, all this in a local definition of the _energy_changed() method. The example could even include a feature that enables/disables the self._calc.energy = value (and subsequent update) code based on instrument's choice (probably yet another PV).

[prjemian]

Per our discussion in the weekly coalition teleconferences, let's use the pint package to manage the units conversion. Example can show both methods (1 and 2) from above.

[ambarb]

should 2 also include hand set energy offset as an attribute that feeds into hklpy_object.calc.energy to take care of https://github.com/bluesky/hklpy/issues/25? Beamlines can decide if it should point to a PV or not.

[prjemian]

yes

[prjemian]

simple using pint for units:

(bluesky_2020_9) prjemian@poof ~/.../projects/apstools $ ipython
Python 3.8.2 (default, Mar 26 2020, 15:53:00) 
Type 'copyright', 'credits' or 'license' for more information
IPython 7.16.1 -- An enhanced Interactive Python. Type '?' for help.

In [1]: import pint                                                                                                                                                                                      

In [2]: ureg = pint.UnitRegistry()                                                                                                                                                                       

In [3]: value = 931 * ureg.eV                                                                                                                                                                            

In [4]: value.to(ureg.keV)                                                                                                                                                                               
Out[4]: 0.931 <Unit('kiloelectron_volt')>

In [5]: print((931 * ureg.eV).to(ureg.keV))                                                                                                                                                              
0.931 kiloelectron_volt

In [6]: value._magnitude                                                                                                                                                                                 
Out[6]: 931

In [7]:                                                                                                                                                                                                  

BUT, convert eV to angstrom, does not work (pint does not know the physics):

In [7]: value.to(ureg.angstrom)                                                                                                                                                                          
---------------------------------------------------------------------------
DimensionalityError                       Traceback (most recent call last)
<ipython-input-7-136c6af86708> in <module>
----> 1 value.to(ureg.angstrom)

~/Apps/anaconda/envs/bluesky_2020_9/lib/python3.8/site-packages/pint/quantity.py in to(self, other, *contexts, **ctx_kwargs)
    616         other = to_units_container(other, self._REGISTRY)
    617 
--> 618         magnitude = self._convert_magnitude_not_inplace(other, *contexts, **ctx_kwargs)
    619 
    620         return self.__class__(magnitude, other)

~/Apps/anaconda/envs/bluesky_2020_9/lib/python3.8/site-packages/pint/quantity.py in _convert_magnitude_not_inplace(self, other, *contexts, **ctx_kwargs)
    565                 return self._REGISTRY.convert(self._magnitude, self._units, other)
    566 
--> 567         return self._REGISTRY.convert(self._magnitude, self._units, other)
    568 
    569     def _convert_magnitude(self, other, *contexts, **ctx_kwargs):

~/Apps/anaconda/envs/bluesky_2020_9/lib/python3.8/site-packages/pint/registry.py in convert(self, value, src, dst, inplace)
    942             return value
    943 
--> 944         return self._convert(value, src, dst, inplace)
    945 
    946     def _convert(self, value, src, dst, inplace=False, check_dimensionality=True):

~/Apps/anaconda/envs/bluesky_2020_9/lib/python3.8/site-packages/pint/registry.py in _convert(self, value, src, dst, inplace)
   1802                 value, src = src._magnitude, src._units
   1803 
-> 1804         return super()._convert(value, src, dst, inplace)
   1805 
   1806     def _get_compatible_units(self, input_units, group_or_system):

~/Apps/anaconda/envs/bluesky_2020_9/lib/python3.8/site-packages/pint/registry.py in _convert(self, value, src, dst, inplace)
   1408 
   1409         if not (src_offset_unit or dst_offset_unit):
-> 1410             return super()._convert(value, src, dst, inplace)
   1411 
   1412         src_dim = self._get_dimensionality(src)

~/Apps/anaconda/envs/bluesky_2020_9/lib/python3.8/site-packages/pint/registry.py in _convert(self, value, src, dst, inplace, check_dimensionality)
    975             # then the conversion cannot be performed.
    976             if src_dim != dst_dim:
--> 977                 raise DimensionalityError(src, dst, src_dim, dst_dim)
    978 
    979         # Here src and dst have only multiplicative units left. Thus we can

DimensionalityError: Cannot convert from 'electron_volt' ([length] ** 2 * [mass] / [time] ** 2) to 'angstrom' ([length])

[gfabbris]

Keeping in mind that I haven't used hklpy at 4IDD yet (only a few times at CSX), here is my take:

• It'd be good to have the diffractometer device to report both energy and wavelength, because this turns out to be a bit of a personal preference. Generally, people with a scattering background prefer to work with wavelength, but those with a spectroscopy background prefer energy.

• As I understand, the energy calibration issue occurs when there are energy drifts in the middle of an experiment, correct? This is typically not an issue at 4IDD, we just calibrate the energy in the beginning of the experiment (by resetting the mono angle). But I suggest that any energy offsets should stay within the monochromator device, and only pass the "final" energy to hklpy.

• While method 1 seems to be better to me as a standard operation, it'd be good to have the option to manually set the energy used by hklpy, because sometimes we want to calculate positions using a different energy.

[prjemian]

Here's an example that an instrument could use as a starting point:

import gi
gi.require_version('Hkl', '5.0')
import hkl.diffract
import pint

ureg = pint.UnitRegistry()

class LocalDiffractometer(hkl.diffract.Diffractometer):
    # ...
    energy = Component(EpicsSignal, "EPICS:diffractometer:energy")
    energy_EGU = Component(EpicsSignal, "EPICS:diffractometer:energy.EGU")
    energy_update_calc = Component(EpicsSignal, "EPICS:diffractometer:energy_update_hkl")
    # ...

    def _energy_changed(self, value=None, **kwargs):
        '''
        Callback indicating that the energy signal was updated
        '''
        if energy_update_calc.get() in (1, "Yes", "locked", "OK"):
            units = self.energy_EGU.get()
            logger.debug(
                '{0.name} energy changed: {1} {2}'.format(
                    self, value, units))
            keV = pint.Quantity(value, units).to(ureg.keV)
            self._calc.energy = keV
            self._update_position()

[prjemian]

If 2, where the instrument always uses the same units for energy (pretty typical or the beamline controls would have other problems), then do not define self.energy_EGU and use like this instead (if conversion is even needed):

        units = "eV"

[prjemian]

@gfabbris :

It'd be good to have the diffractometer device to report both energy and wavelength

Are you referring to the wh() report? (This is from the apstools.diffractometer.DiffractometerMixin() class.)

[prjemian]

The point of leaving out any mention of energy from the hklpy code is to allow it to be agnostic to the probe type (primarily: X-ray vs neutron). BUT, neutron users could always override these methods in hkl.calc.CalcRecip():

    @property
    def wavelength(self):
        '''The wavelength associated with the geometry, in angstrom'''
        return self._geometry.wavelength_get(self._units)

    @wavelength.setter
    def wavelength(self, wavelength):
        self._geometry.wavelength_set(wavelength, self._units)

    @property
    def energy(self):
        '''The energy associated with the geometry, in keV'''
        return A_KEV / self.wavelength

    @energy.setter
    def energy(self, energy):
        self.wavelength = A_KEV / energy

[ambarb]

I do not think we should leave energy out of hklpy for reasons described in point 1 response. One should be able to chose if you want to use energy or wavelength, and the documentation should be clear on this.

To reply to @gfabbris 3 bulltet points: point 1 It seems to be both configuration attributes of hklpy object are possible to include. The beamline typically chooses what to record in the documents sent to the database. Claudio and Stuart made these decisions, and it is energy because this is how we control the incident beamline energy, not via wavelength. For REXS, it is probably nicer to use energy and not wavelength because maybe you want to tune to a few different energies that have specific features for a given atom. Therefore, I think it is important to have both options and not make these options "technique specific" (x-rays versus neturons versus electrons, diffraction versus spectroscopy version spectroscopic diffraction). The neutron instruments I have used have "banks" of energy or setting the incident wavelength at very specific options. If finer wavelength control is possible, it is so complex there is an entire system that the beamline scientists and controls have set up to deliver fully processed data. But perhaps we can get @stuartcampbell to weigh in here.

point 2 What do you mean by monochrometer device? The ophyd object that is the monochrometer? For CSX, the energy resolution is not high, so we do not see energy drifts as a function of "time" but more so as a function of accelerator configuration, switching between gratings, point of energy focus (which is configurable at CSX), and other complex factors that we are trying to remove. But in this sense, the offset is NOT stable in time.

If monochrometer = mono IOC, then I think this is too big of a requirement to force on people. If you mean in the ophyd device, then this affects several beamlines and then this requirement for hklpy compatible monochrometer devices needs to be very clear, and we need to provide an example in documentation how to do this. My only concern is switching between energies. The energy calibration is not constant or linear (say O Kedge at 530eV versus Cu L3 at 930eV). This could create a huge mess if these things are not managed carefully. If we decide to manage in hklpy, I think the impact of this "mess" could be less and easier to sort out.

point 3 Method 1 sounds nice and very convenient, and was my first choice. However, Method 1 also requires more complex monochrometer ophyd devices and a more complex entry of nm verses angstroms for unit cell parameters and all of the metadata associated with this. Now, people are discussing making sample databases require this information so it can be automatically added to scan metadata. I think coding up Method 1 to "dynamically" switch between whatever a user wants to do would be more complex, where as method 2 is a configuration decided by beamline scientist to work best with the workflows they setup for metadata pre-scan, plans needed, and automated analysis. For practical purposes, I think Method 2 represents something we can do now that would lay the ground work for Method 1.

[ambarb]

I also feel like this discussion is affecting larger groups of people and that we need to come up with requirements documents from the APS and NSLS-II and anyone else who cares. And then work from those.

If we agree, then I can start getting people's requirement at NSLS-II.

[prjemian]

Good idea. If you know the GitHub user name of interested parties, invite them to contribute to the discussion here. Such as: @YChoi

[ambarb]

most of them don't use git (or just ignore git emails). I think the default setting is actually to have notifications of [at]so-and-so is OFF

[prjemian]

@junspark @Jasonzhouquan @cpchuang @rodolakis You are invited to comment on how energy or wavelength are used in your use of diffractometer (hkl) operations. This discussion centers on how to interface the beam line's sens of energy with the wavelength used by the calculation converting between (hkl) and motor positions.

Invite others to participate as you wish (send emails if they do not use github yet).

[ambarb]

another comment to @gfabbris point 2: I've used beamlines at the APS that don't know the beamline energy so exactly that diffraction just works. To deal with that, I have done all of the following:

• using reflex to fudge a b c alpha beta gamma once you have found a few Bragg peaks and you don't have any phase transitions that rely on knowing how a compares to c (for example)
• don't use whole numbers for H, K, L and just accept H=1.012 (or whatever it is)
• adjust the lattice constants by hand to make the orientation matrix work

If I had a way to determine the energy in those cases, I would have preferred to do that. In @gfabbris point, we should leave people options to do the experiment the way they want to do the experiment (and at beamline scientists discretion in how things work).

[gfabbris]

Jut to clarify, in point 2 I meant the monochromator ophyd object. But in any case, I agree with @ambarb point that it may be better to leave the mono settings alone if you need to change the offset rather often, and as far as I can tell it'd be rather straightforward to have an offset in hklpy that can be just ignored if it's not needed.

[prjemian]

an additive energy_offset (the units are specific to the instrument) such as:

    energy_offset = Signal(value=0)

then

    def _energy_changed(self, value=None, **kwargs):
        '''
        Callback indicating that the energy signal was updated
        '''
        if energy_update_calc.get() in (1, "Yes", "locked", "OK"):
            local_energy = value + self.energy_offset.get()
            units = self.energy_EGU.get()   # or fixed, such as:  units = "eV"
            logger.debug(
                '{0.name} energy changed: {1} {2}'.format(
                    self, value, units))
            keV = pint.Quantity(local_energy, units).to("keV")
            self._calc.energy = keV
            self._update_position()

[ambarb]

@prjemian Are we only collecting requirements in relationship to the proposed new class? Or are we collecting requirements for all things. Or new class now, and all other things later?

If class now, do you want to fold in the other energy related Issues?

• https://github.com/bluesky/hklpy/issues/30
• https://github.com/bluesky/hklpy/issues/26
• https://github.com/bluesky/hklpy/issues/25

[prjemian]

Good question. This issue was started specific to wavelength which seems to be an unpopular suggestion. A way to resolve #26 (for energy) was just shown. Should we keep the issue open but change it's title as you suggest? We could open issues gathered in due course.

[prjemian]

Per suggestion above, see https://github.com/BCDA-APS/apstools/pull/433.

[JZTischler]

Diffraction, of course, uses wavelength, but x-ray people usually like to use energy. The BlueSky end should be able to use either and have an internal routine for converting back and forth wavelength<--> energy (if you set one the other is calculated). This routine will have to be customized for photons, electron, and neutrons. To set or read from a physical monochromator, just pick either wavelength or energy as you prefer, the monochromator needs to have the same conversion routines. Since hc is now exact, and the electron and neutron masses are know to a part in 1e-10, there should be no difference between the BlueSkY conversion and the monochromator's conversion.

The only other think that I would add, is that it would be nice to be able to establish multiple matrices. E.g. one matrix for the substrate, and another for the film. Or one for each of the 4 domains that might form during a phase transformation.

[prjemian]

@JZTischler : Thanks for the input. This code keeps a list of samples (in the background), each with their own Lattice definition and orientation matrix. One of those is selected to be the primary at any time. It is possible to change between them by the text name they were given when first defined.

[prjemian]

Summary

• Will not add wavelength to Diffractometer class.
• Energy with units and offset is supported in Diffractometer class.
• Monochromator energy can be connected to Diffractometer energy, as shown above.
• Conversion between X-ray energy and wavelength is performed in methods of the calc module, as shown above.
• For use with neutrons (or other radiation sources), the Diffractometer class can be subclassed and the conversion methods replaced as appropriate.
• Energy is reported in both the diffract.Diffractometer.pa() and diffract.Diffractometer.wh() reports.

Thanks to all who contributed to this discussion.