From d18a0311a80cb57a90f21a7bf4a46ee944e05369 Mon Sep 17 00:00:00 2001
From: Duc Le <duc.le@stfc.ac.uk>
Date: Sat, 26 May 2018 00:32:35 +0100
Subject: [PATCH] Re #22359 - changed code to use an input yaml instrument file
Refactored code ISISDisk.py and ISISFermi.py into Instruments.py
Deprecated ISISDisk and ISISFermi
Changed GUI to use Instruments
---
scripts/PyChop/ISISDisk.py | 3 +-
scripts/PyChop/ISISFermi.py | 2 +-
scripts/PyChop/Instruments.py | 798 ++++++++++++++++++++++++++++++++++
scripts/PyChop/MulpyRep.py | 25 +-
scripts/PyChop/PyChop2.py | 3 +
scripts/PyChop/PyChopGui.py | 22 +-
scripts/PyChop/__init__.py | 2 +-
scripts/PyChop/let.yaml | 138 ++++++
scripts/PyChop/maps.yaml | 130 ++++++
9 files changed, 1100 insertions(+), 23 deletions(-)
create mode 100644 scripts/PyChop/Instruments.py
create mode 100644 scripts/PyChop/let.yaml
create mode 100644 scripts/PyChop/maps.yaml
diff --git a/scripts/PyChop/ISISDisk.py b/scripts/PyChop/ISISDisk.py
index c9d0bfdf4a6..c73c86acc3f 100644
--- a/scripts/PyChop/ISISDisk.py
+++ b/scripts/PyChop/ISISDisk.py
@@ -20,7 +20,7 @@ class ISISDisk:
def __init__(self, instname=None, variant=None, freq=None):
warnings.warn("The ISISDisk class is deprecated and will be removed in the next Mantid version. "
- "Please use the Instrument class or the official PyChop2 CLI interface.", DeprecationWarning)
+ "Please use the Instrument class or the official PyChop CLI interface.", DeprecationWarning)
if instname:
self.setInstrument(instname, variant)
self.freq = 0
@@ -465,6 +465,7 @@ class ISISDisk:
dist, samDist, DetDist, fracEi = tuple([self.dist, self.chop_samp, self.samp_det, self.frac_ei])
modSamDist = dist[-1] + samDist
totDist = modSamDist + DetDist
+ print(chop_times)
for i in range(len(dist)):
plt.plot([-20000, 120000], [dist[i], dist[i]], c='k', linewidth=1.)
for j in range(len(chop_times[i][:])):
diff --git a/scripts/PyChop/ISISFermi.py b/scripts/PyChop/ISISFermi.py
index 950eef7b8f6..1b49160c4ab 100644
--- a/scripts/PyChop/ISISFermi.py
+++ b/scripts/PyChop/ISISFermi.py
@@ -85,7 +85,7 @@ class ISISFermi:
def __init__(self, instname=None, choppername='', freq=0):
warnings.warn("The ISISFermi class is deprecated and will be removed in the next Mantid version. "
- "Please use the Instrument class or the official PyChop2 CLI interface.", DeprecationWarning)
+ "Please use the Instrument class or the official PyChop CLI interface.", DeprecationWarning)
if instname:
self.setInstrument(instname, choppername, freq)
self.diskchopper_phase = self.__DiskChopperMode[instname]
diff --git a/scripts/PyChop/Instruments.py b/scripts/PyChop/Instruments.py
new file mode 100644
index 00000000000..975eccae0d3
--- /dev/null
+++ b/scripts/PyChop/Instruments.py
@@ -0,0 +1,798 @@
+"""
+This module is a wrapper around a set of instrument parameters (to be read from a YAML file)
+and methods which then call either Chop.py or MulpyRep.py to do the resolution calculations.
+"""
+
+from __future__ import (absolute_import, division, print_function)
+from six import string_types
+import numpy as np
+import yaml
+import warnings
+import copy
+try:
+ from . import Chop, MulpyRep
+except:
+ import Chop
+ import MulpyRep
+from scipy.interpolate import interp1d
+from scipy.special import erf
+from scipy import constants
+
+# Some global constants
+SIGMA2FWHM = 2 * np.sqrt(2 * np.log(2))
+SIGMA2FWHMSQ = SIGMA2FWHM ** 2
+E2V = np.sqrt((constants.e / 1000) * 2 / constants.neutron_mass) # v = E2V * sqrt(E) veloc in m/s, E in meV
+E2L = 1.e23 * constants.h**2 / (2 * constants.m_n * constants.e) # lam = sqrt(E2L / E) lam in Angst, E in meV
+E2K = constants.e * 2 * constants.m_n / constants.hbar**2 / 1e23 # k = sqrt(E2K * E) k in 1/Angst, E in meV
+
+def wrap_attributes(obj, inval, allowed_var_names):
+ for key in allowed_var_names:
+ if inval:
+ if hasattr(inval, key):
+ setattr(obj, key, getattr(inval, key))
+ elif hasattr(inval, '__getitem__') and key in inval:
+ setattr(obj, key, inval[key])
+
+def argparser(args, kwargs, argnames, defaults=None):
+ argdict = {key: val for key, val in zip(argnames, defaults if defaults else [None] * len(argnames))}
+ for key in kwargs:
+ argdict[key] = kwargs[key]
+ for idx in range(len(args)):
+ argdict[argnames[idx]] = args[idx]
+ return argdict
+
+def soft_hat(x, p):
+ """
+ ! Soft hat function, from Herbert subroutine library.
+ ! For rescaling t-mod at low energy to account for broader moderator term
+ """
+ x = np.array(x)
+ sig2fwhh = np.sqrt(8*np.log(2))
+ height, grad, x1, x2 = tuple(p[:4])
+ sig1, sig2 = tuple(np.abs(p[4:6]/sig2fwhh))
+ # linearly interpolate sig for x1<x<x2
+ sig = ((x2-x)*sig1-(x1-x)*sig2)/(x2-x1)
+ if np.shape(sig):
+ sig[x < x1] = sig1
+ sig[x > x2] = sig2
+ # calculate blurred hat function with gradient
+ e1 = (x1-x) / (np.sqrt(2)*sig)
+ e2 = (x2-x) / (np.sqrt(2)*sig)
+ y = (erf(e2)-erf(e1)) * ((height+grad*(x-(x2+x1)/2))/2)
+ y = y + 1
+ return y
+
+
+class FermiChopper(object):
+ """
+ Class which represents a Fermi chopper package
+ """
+
+ __allowed_var_names = ['name', 'pslit', 'pslat', 'radius', 'rho', 'tjit', 'fluxcorr']
+
+ def __init__(self, inval=None):
+ wrap_attributes(self, inval, self.__allowed_var_names)
+
+ def __repr__(self):
+ return self.name if self.name else 'Undefined Fermi chopper package'
+
+ def getWidthSquared(self, Ei, freq):
+ return Chop.tchop(freq, Ei, self.pslit / 1000., self.radius / 1000., self.rho / 1000.)
+
+ def getWidth(self, *args):
+ """ Calculates the chopper time width in seconds for a given neutron energy (Ei) """
+ return np.sqrt(self.getWidthSquared(*args))
+
+ def getTransmission(self, Ei, freq):
+ """ Calculates the chopper transmission """
+ dslat = (self.pslit + self.pslat) / 1000
+ return Chop.achop(Ei, freq, dslat, self.pslit / 1000., self.radius / 1000., self.rho / 1000.) / dslat
+
+
+class ChopperSystem(object):
+ """
+ Class which represents a set (list) of choppers in a line
+ """
+
+ __allowed_var_names = ['name', 'chop_sam', 'sam_det', 'aperture_width', 'aperture_height', 'choppers', 'variants',
+ 'frequency_matrix', 'constant_frequencies', 'max_frequencies', 'default_frequencies',
+ 'source_rep', 'emission_time', 'overlap_ei_frac', 'ei_limits', 'flux_ref_slot', 'flux_ref_freq']
+
+ def __init__(self, inval=None):
+ # Default values
+ self.source_rep = 50
+ self.emission_time = 0
+ self.overlap_ei_frac = 0.9
+ self._ei = 22
+ # Parse input values (if any)
+ wrap_attributes(self, inval, self.__allowed_var_names)
+ self._parse_choppers()
+ self._parse_variants()
+ self.phase = self.defaultPhase
+ self.frequency = self.default_frequencies
+
+ def __repr__(self):
+ return self.name if self.name else 'Undefined disk chopper system'
+
+ def _parse_choppers(self):
+ """Parses the choppers list to determine how to handle resolution and flux calculations"""
+ if not self.choppers:
+ return
+ self.distance = []
+ self.nslot = []
+ self.slot_ang_pos = []
+ self.slot_width = []
+ self.guide_width = []
+ self.radius = []
+ self.numDisk = []
+ self.isFermi = False
+ self.isPhaseIndependent = []
+ self.defaultPhase = []
+ for idx, chopper in enumerate(self.choppers):
+ self.distance.append(chopper['distance'])
+ if 'packages' in chopper:
+ self.isFermi = True
+ self.idxFermi = idx
+ self.packages = {key: FermiChopper(val) for key, val in list(chopper['packages'].items())}
+ self.nslot.append(1) # Assume Fermi chopper is curved, will not transmit PI pulse.
+ self.slot_ang_pos.append(None)
+ self.slot_width.append(10.)
+ self.guide_width.append(10.)
+ self.radius.append(290.)
+ self.numDisk.append(1)
+ else:
+ self.nslot.append(
+ chopper['nslot'] if ('nslot' in chopper and chopper['nslot']) else len(chopper['slot_ang_pos']))
+ self.slot_ang_pos.append(chopper['slot_ang_pos'] if ('slot_ang_pos' in chopper) else None)
+ self.slot_width.append(chopper['slot_width'])
+ self.guide_width.append(chopper['guide_width'])
+ self.radius.append(chopper['radius'])
+ self.numDisk.append(2 if ('isDouble' in chopper and chopper['isDouble']) else 1)
+ self.isPhaseIndependent.append(
+ True if ('isPhaseIndependent' in chopper and chopper['isPhaseIndependent']) else False)
+ self.defaultPhase.append(chopper['defaultPhase'] if 'defaultPhase' in chopper else None)
+ if not any(self.slot_ang_pos):
+ self.slot_ang_pos = None
+ self.isPhaseIndependent = [ii for ii in range(len(self.isPhaseIndependent)) if self.isPhaseIndependent[ii]]
+ self.defaultPhase = [self.defaultPhase[ii] for ii in self.isPhaseIndependent]
+ self._instpar = [self.distance, self.nslot, self.slot_ang_pos, self.slot_width, self.guide_width, self.radius,
+ self.numDisk, self.sam_det, self.chop_sam, self.source_rep, self.emission_time,
+ self.overlap_ei_frac, self.isPhaseIndependent]
+ if self.isFermi:
+ self.package = self.packages.keys()[0]
+
+ def _parse_variants(self):
+ if 'variants' not in self.__dict__:
+ return
+ variant_keys = []
+ for var in self.variants:
+ if ('default' in self.variants[var].keys() and self.variants[var]['default']) or var is None:
+ self._default_variant = self._variant = var
+ if var:
+ [variant_keys.append(key) for key in self.variants[var].keys() if 'default' not in key]
+ self._variant_defaults = {}
+ for key in set(variant_keys):
+ self._variant_defaults[key] = copy.deepcopy(getattr(self, key))
+ if '_variant' not in self.__dict__:
+ self._default_variant = self.variants.keys()[0]
+ warnings.warn('No default variants defined. Using ''%s'' as default' % (self._default_variant), SyntaxWarning)
+ self.variant = self._default_variant
+
+ # Define getters/setters here to be backwards compatible with old PyChop2. Actually use properties underneath
+ def setChopper(self, *args, **kwargs):
+ """
+ Set the chopper package type (Fermi instruments) or variant (LET).
+
+ maps = Instrument('MAPS')
+ maps.setChopper('A', 400) # Sets package A at 400 Hz.
+ maps.setChopper(package='A', freq=400) # Explicit keywords
+ let = Instrument('LET')
+ let.setChopper('High resolution', [280, 140]) # Change to the high resolution variant at 280 Hz
+ let.setChopper(variant='High resolution')
+ """
+ argdict = argparser(args, kwargs, ['package' if self.isFermi else 'variant', 'freq'])
+ if self.isFermi:
+ self.package = argdict['package']
+ if hasattr(self, 'variants') and argdict['package'] in self.variants:
+ self.variant = argdict['package']
+ elif argdict['variant']:
+ self.variant = argdict['variant']
+ if argdict['freq']:
+ self.frequency = argdict['freq']
+
+ def getChopper(self):
+ return self.package if self.isFermi else self.variant
+
+ def setFrequency(self, *args, **kwargs):
+ """
+ Set the chopper frequency(ies) and (optionally) phase(s).
+
+ maps = Instrument('MAPS', 'A')
+ maps.setFrequency(400) # Sets frequency to 400 Hz.
+ maps.setFrequency([400, 100], 1) # Sets Fermi to 400Hz, disk to 100Hz, and multi-rep mode
+ maps.setFrequency(freq=[400, 100], phase=1)
+ let = Instrument('LET')
+ let.setFrequency([240, 120]) # Sets resolution chopper to 240Hz, pulse removal to 120Hz
+ let.setFrequency([240, 120], phase=-20000) # Additionally sets the frame overlap phase to -20000us
+ """
+ argdict = argparser(args, kwargs, ['freq', 'phase'])
+ if argdict['freq']:
+ self.frequency = argdict['freq']
+ if argdict['phase']:
+ self.phase = argdict['phase']
+
+ def getFrequency(self):
+ return self.frequency
+
+ def setEi(self, Ei):
+ """Sets the (focussed) incident energy"""
+ self.ei = Ei
+
+ def getEi(self):
+ return self.ei
+
+ def getAllowedEi(self, Ei_in=None):
+ return self._MulpyRepDriver(Ei_in, calc_res=False)[0]
+
+ def plotMultiRepFrame(self, h_plt=None, Ei_in=None, frequency=None):
+ """
+ Plots the time-distance diagram into a given Matplotlib axes, i
+ for a give focused incident energy (in meV) and chopper frequencies (in Hz).
+ """
+ if h_plt is None:
+ try:
+ from matplotlib import pyplot
+ except ImportError:
+ raise RuntimeError('plotMultiRepFrame: Cannot import matplotlib')
+ plt = pyplot
+ else:
+ plt = h_plt
+ Ei = self.ei if Ei_in is None else Ei_in
+ if Ei is None:
+ raise ValueError('Incident energy has not been specified')
+ if frequency:
+ oldfreq = self.freq
+ self.setFrequency(frequency)
+ Eis, _, _, lines, chop_times = tuple(self._MulpyRepDriver(Ei_in, calc_res=False))
+ if frequency:
+ self.setFrequency(oldfreq)
+ modSamDist = self.distance[-1] + self.chop_sam
+ totDist = modSamDist + self.sam_det
+ for i in range(len(self.distance)):
+ plt.plot([-20000, 120000], [self.distance[i], self.distance[i]], c='k', linewidth=1.)
+ for j in range(len(chop_times[i][:])):
+ plt.plot(chop_times[i][j], [self.distance[i], self.distance[i]], c='white', linewidth=1.)
+ plt.plot([-20000, 120000], [totDist, totDist], c='k', linewidth=2.)
+ for i in range(len(lines)):
+ x0 = (-lines[i][0][1] / lines[i][0][0] - lines[i][1][1] / lines[i][1][0]) / 2.
+ x1 = ((modSamDist-lines[i][0][1]) / lines[i][0][0] + (modSamDist-lines[i][1][1]) / lines[i][1][0]) / 2.
+ plt.plot([x0, x1], [0, modSamDist], c='b')
+ x2 = ((totDist-lines[i][0][1]) / lines[i][0][0] + (totDist-lines[i][1][1]) / lines[i][1][0]) / 2.
+ lineM = totDist / x2
+ plt.plot([x1, x2], [modSamDist, totDist], c='b')
+ newline = [lineM * np.sqrt(1 + self.overlap_ei_frac), modSamDist - lineM * np.sqrt(1 + self.overlap_ei_frac) * x1]
+ x3 = (totDist-newline[1]) / (newline[0])
+ plt.plot([x1, x3], [modSamDist, totDist], c='r')
+ newline = [lineM * np.sqrt(1 - self.overlap_ei_frac), modSamDist - lineM * np.sqrt(1 - self.overlap_ei_frac) * x1]
+ x4 = (totDist-newline[1]) / (newline[0])
+ plt.plot([x1, x4], [modSamDist, totDist], c='r')
+ plt.text(x2, totDist+0.2, "{:3.1f}".format(Eis[i]))
+ xmax = (1.e6 / self.source_rep) if hasattr(self, 'source_rep') and self.source_rep else np.max(chop_times)
+ if h_plt is None:
+ plt.xlim(0, xmax)
+ plt.xlabel(r'TOF ($\mu$sec)')
+ plt.ylabel('Distance (m)')
+ plt.show()
+ else:
+ plt.set_xlim(0, xmax)
+ plt.set_xlabel(r'TOF ($\mu$sec)')
+ plt.set_ylabel(r'Distance (m)')
+
+ def getWidthSquared(self, Ei_in=None):
+ return self.getWidth(Ei_in, squared=True)
+
+ def getWidth(self, Ei_in=None, squared=False):
+ """Returns the chopper time width (FWHM) at the (final) chopper in microseconds"""
+ Ei = self.ei if Ei_in is None else Ei_in
+ if self.isFermi:
+ return self._ChopDriver(Ei_in, squared), None
+ else:
+ chop_times = self._MulpyRepDriver(Ei_in, calc_res=False)[1]
+ # Output of MulpyRep is FWHM in us - want it in seconds for later calculations
+ wd = ((chop_times[1][1] - chop_times[1][0]) / 2. / 1.e6, (chop_times[0][1] - chop_times[0][0]) / 2. / 1.e6)
+ return (wd[0]**2, wd[1]**2) if squared else wd
+
+ def getDistances(self):
+ """ Returns the (mod->final_chop, aperture->final, chop->sam, sam->det, mod-first_chop) distances for instrument """
+ mod_chop = self.choppers[-1]['distance']
+ ap_chop = self.choppers[-1]['aperture_distance'] if ('aperture_distance' in self.choppers[-1]) else mod_chop
+ return (mod_chop, ap_chop, self.chop_sam, self.sam_det, self.choppers[0]['distance'])
+
+ def getTransmission(self, Ei_in=None, frequency=None, hires=False):
+ """ Calculates the flux transmission fraction through the chopper system at specified Ei and frequency """
+ Ei = self.ei if Ei_in is None else Ei_in
+ freq = frequency if frequency is not None else self._long_frequency[-1]
+ if self.isFermi:
+ x0, x1 = (self.choppers[-1]['distance'], self.chop_sam)
+ magic = (84403.06 / x0 / (x1 + x0)) # Some magical conversion factor (??)
+ fudge = self.packages[self.package].fluxcorr # A chopper package dependent fudge factor
+ return self.packages[self.package].getTransmission(Ei, freq) * magic / fudge
+ else:
+ return (self.slot_width[-1] / self.flux_ref_slot) * (self.flux_ref_freq / freq)
+
+ def _get_state(self, Ei_in=None):
+ return hash((self.variant, self.package, tuple(self.frequency), tuple(self.phase), Ei_in if Ei_in else self.ei))
+
+ def _removeLowIntensityReps(self, Eis, lines, Ei=None):
+ # Removes reps with Ei where there are no neutrons
+ if not hasattr(self, 'ei_limits') or not self.ei_limits:
+ return Eis, lines
+ Eis = np.array(Eis)
+ idx = (Eis >= self.ei_limits[0]) * (Eis <= self.ei_limits[1])
+ # Always keeps desired rep even if outside of range
+ if Ei:
+ idx += ((np.abs(Eis - Ei) / np.abs(Eis)) < 0.1)
+ Eis = Eis[idx]
+ lines = np.array(lines)[idx]
+ return Eis, lines
+
+ def _MulpyRepDriver(self, Ei_in=None, calc_res=True):
+ """Private method to calculate resolution for given Ei from chopper opening times"""
+ Ei = self.ei if Ei_in is None else Ei_in
+ if '_saved_state' not in self.__dict__ or (self._saved_state[0] != self._get_state(Ei)):
+ Eis, all_times, chop_times, lastChopDist, lines = MulpyRep.calcChopTimes(Ei, self._long_frequency, self._instpar, self.phase)
+ Eis, lines = self._removeLowIntensityReps(Eis, lines, Ei)
+ self._saved_state = [self._get_state(Ei), Eis, chop_times, lastChopDist, lines, all_times]
+ else:
+ Eis, chop_times, lastChopDist, lines, all_times = tuple(self._saved_state[1:])
+ if calc_res:
+ res_el, percent, chop_width, mod_width = MulpyRep.calcRes(Eis, chop_times, lastChopDist, self.chop_sam, self.sam_det)
+ return res_el, percent, chop_width, mod_width
+ else:
+ return [Eis, chop_times, lastChopDist, lines, all_times]
+
+ def _ChopDriver(self, Ei_in=None, squared=False):
+ """Private method to calculate resolution for given Ei from Fermi chopper"""
+ Ei = self.ei if Ei_in is None else Ei_in
+ if squared:
+ return self.packages[self.package].getWidthSquared(Ei, self._long_frequency[-1]) * SIGMA2FWHMSQ
+ else:
+ return self.packages[self.package].getWidth(Ei, self._long_frequency[-1]) * SIGMA2FWHM
+
+ @property
+ def frequency(self):
+ return self._frequency
+
+ @frequency.setter
+ def frequency(self, value):
+ freq = self.default_frequencies
+ if not hasattr(value, '__len__'):
+ value = [value]
+ freq = [value[i] if i < len(value) else freq[i] for i in range(len(freq))]
+ if self.max_frequencies and not (freq <= self.max_frequencies):
+ raise ValueError('Value of frequencies outside maximum allowed')
+ self._frequency = freq
+ if hasattr(self, 'constant_frequencies') and self.constant_frequencies:
+ f0 = self.constant_frequencies
+ else:
+ f0 = [0] * np.shape(self.frequency_matrix)[0]
+ self._long_frequency = np.dot(self.frequency_matrix, freq) + f0
+
+ @property
+ def phase(self):
+ return self._phase
+
+ @phase.setter
+ def phase(self, value):
+ phase = self.defaultPhase
+ if not hasattr(value, '__len__'):
+ value = [value]
+ self._phase = [value[i] if i < len(value) else phase[i] for i in range(len(phase))]
+
+ @property
+ def ei(self):
+ return self._ei
+
+ @ei.setter
+ def ei(self, value):
+ if value < 0:
+ raise ValueError('Incident neutron energy cannot be less than zero')
+ self._ei = value
+
+ @property
+ def package(self):
+ return self._package if self.isFermi else None
+
+ @package.setter
+ def package(self, value):
+ if not self.isFermi:
+ raise AttributeError('Cannot set Fermi chopper package on this instrument')
+ if value not in self.packages.keys():
+ raise ValueError('Fermi package ''%s'' not recognised. Allowed values are: %s'
+ % (value, ', '.join(self.packages.keys())))
+ self._package = value
+
+ @property
+ def variant(self):
+ return self._variant if hasattr(self, 'variants') and self.variants else None
+
+ def getAllowedChopper(self):
+ return self.packages.keys() if self.isFermi else (self.variants.keys() if self.variants else None)
+
+ @variant.setter
+ def variant(self, value):
+ if 'variants' not in self.__dict__:
+ raise AttributeError('This instrument has no variants to set')
+ for prop in self._variant_defaults:
+ setattr(self, prop, copy.deepcopy(self._variant_defaults[prop]))
+ self._variant = value
+ if value not in self.variants.keys():
+ raise ValueError('Variant ''%s'' not recognised. Allowed values are: %s' % (value, ', '.join(self.variants.keys())))
+ for prop in self.variants[value]:
+ if prop == 'choppers':
+ for idx, chopper in enumerate(self.variants[value][prop]):
+ if chopper:
+ for ky in chopper:
+ self.choppers[idx][ky] = chopper[ky]
+ elif prop in self.__allowed_var_names:
+ setattr(self, prop, self.variants[value][prop])
+ self._parse_choppers()
+
+ @property
+ def tjit(self):
+ return self.packages[self.package].tjit if self.isFermi else 0.
+
+ @tjit.setter
+ def tjit(self, value):
+ self.packages[self.package].tjit = value
+
+
+class Moderator(object):
+ """
+ Class which represents a neutron moderator
+ """
+
+ __allowed_var_names = ['name', 'imod', 'ch_mod', 'mod_pars', 'mod_scale_fn', 'mod_scale_par', 'theta',
+ 'source_rep', 'emission_time', 'measured_flux']
+
+ def __init__(self, inval=None):
+ wrap_attributes(self, inval, self.__allowed_var_names)
+ if hasattr(self, 'measured_flux') and 'scale_factor' in self.measured_flux:
+ self.measured_flux['flux'] = np.array(self.measured_flux['flux']) * float(self.measured_flux['scale_factor'])
+
+ def __repr__(self):
+ return self.name if self.name else 'Undefined neutron moderator'
+
+ def getWidthSquared(self, Ei):
+ """ Returns the squared time gaussian sigma width due to the sample in s^2 """
+ if self.imod == 0:
+ # CHOP outputs the Gaussian sigma^2 in s^2, we want FWHM^2 in s^2
+ tsqmod = Chop.tchi(self.mod_pars / 1000., Ei) * SIGMA2FWHMSQ
+ elif self.imod == 1:
+ tsqmod = Chop.tikeda(tuple(self.mod_pars), Ei) * SIGMA2FWHMSQ
+ elif self.imod == 2:
+ d0 = self.mod_pars[0]
+ if hasattr(self, 'mod_scale_fn') and self.mod_scale_fn:
+ try:
+ d0 *= globals()[self.mod_scale_fn](Ei, self.mod_scale_par)
+ except KeyError:
+ pass
+ tsqmod = Chop.tchi_2(d0 / 1000., self.mod_pars[1] / 1000., Ei) * SIGMA2FWHMSQ
+ elif self.imod == 3:
+ # Mode for LET - output of polynomial is FWHM in us
+ tsqmod = np.polyval(self.mod_pars, np.sqrt(E2L / Ei))**2 / 1e12
+ else:
+ raise RuntimeError('PyChop: Undefined moderator time profile type %d' % (self.imod))
+ return tsqmod
+
+ def getWidth(self, Ei):
+ """ Calculates the moderator time width in seconds for a given neutron energy (Ei) """
+ if self.imod == 3:
+ # Mode for LET - output of polynomial is FWHM in us
+ return np.polyval(self.mod_pars, np.sqrt(E2L / Ei)) / 1e6
+ else:
+ return np.sqrt(self.getWidthSquared(Ei))
+
+ def getFlux(self, Ei):
+ """ Returns the white beam flux estimate from either measured data (prefered) or analytical model (backup) """
+ return self.getMeasuredFlux(Ei) if (hasattr(self, 'measured_flux') and self.measured_flux) else self.getAnalyticFlux(Ei)
+
+ def getAnalyticFlux(self, Ei):
+ """ Estimate white beam flux from TGP's model of the moderators (ISIS TS1 only) """
+ if all([self.name != modtype for modtype in ['AP', 'CH4', 'H2']]):
+ raise AttributeError('No analytical model for moderator %s' % (self.name))
+ return Chop.flux_calc(np.array(Ei), self.name, self.theta_m * np.pi / 180.)
+
+ def getMeasuredFlux(self, Ei):
+ """ Interpolates flux from a table of measured flux """
+ if not self.measured_flux:
+ raise AttributeError('This instrument does not have a table of measured flux')
+ f = interp1d(self.measured_flux['wavelength'], self.measured_flux['flux'], kind='cubic')
+ return f(np.sqrt(E2L / np.array(Ei)))
+
+ @property
+ def theta_m(self):
+ return self.theta if (hasattr(self, 'theta') and self.theta) else 0.
+
+
+class Sample(object):
+ """
+ Class which represents a sample shape
+ """
+
+ __allowed_var_names = ['name', 'sx', 'sy', 'sz', 'isam', 'gamma']
+
+ def __init__(self, inval=None):
+ wrap_attributes(self, inval, self.__allowed_var_names)
+
+ def __repr__(self):
+ return self.name if self.name else 'Undefined sample'
+
+ def getWidthSquared(self):
+ """ Returns the squared time FWHM due to the sample in s^2 """
+ # At the moment this routine only returns a non-zero y (beam-axis) width
+ return Chop.sam0(self.sx / 1000., self.sy / 1000., self.sz / 1000., self.isam)[1] * SIGMA2FWHMSQ
+
+ def getWidth(self):
+ return np.sqrt(self.getWidthSquared)
+
+ @property
+ def gamma_deg(self):
+ return self.gamma if (hasattr(self, 'gamma') and self.gamma) else 0.
+
+
+class Detector(object):
+ """
+ Class which represents a sample shape
+ """
+
+ __allowed_var_names = ['name', 'idet', 'dd', 'tbin', 'phi']
+
+ def __init__(self, inval=None):
+ wrap_attributes(self, inval, self.__allowed_var_names)
+
+ def __repr__(self):
+ return self.name if self.name else 'Undefined detector'
+
+ def getWidthSquared(self, Ei, en=0):
+ """ Returns the squared time FWHM due to the detector in s^2 """
+ return self.getWidth(Ei, en) ** 2
+
+ def getWidth(self, Ei, en=0):
+ return Chop.detect2(1., 1., np.sqrt(E2K * (Ei-en)), self.idet, self.dd)[3] * SIGMA2FWHM
+
+ @property
+ def phi_deg(self):
+ return self.phi if (hasattr(self, 'phi') and self.phi) else 0.
+
+
+class Instrument(object):
+ """
+ Class which represents a direct geometry neutron spectrometer
+ """
+
+ __allowed_var_names = ['name', 'sample', 'chopper_system', 'moderator', 'detector']
+
+ __child_methods = ['setChopper', 'getChopper', 'setFrequency', 'getFrequency', 'setEi', 'getEi',
+ 'getAllowedEi', 'plotMultiRepFrame']
+
+ __child_properties = ['package', 'variant', 'frequency', 'phase', 'ei', 'tjit']
+
+ def __init__(self, instrument, chopper=None, freq=None):
+ if isinstance(instrument, string_types):
+ try:
+ with open(instrument) as f:
+ instrument = yaml.safe_load(f)
+ except (OSError, IOError) as e:
+ raise RuntimeError('Cannot open file %s . Error is %s' % (instrument, e))
+ if ((hasattr(instrument, 'moderator') or hasattr(instrument, 'chopper_system')) or
+ ('moderator' in instrument or 'chopper_system' in instrument)):
+ wrap_attributes(self, instrument, self.__allowed_var_names)
+ if 'source_rep' in self.moderator:
+ self.chopper_system['source_rep'] = self.moderator['source_rep']
+ if 'emission_time' in self.moderator:
+ self.chopper_system['emission_time'] = self.moderator['emission_time']
+ else:
+ raise RuntimeError('Input to Instrument must be an Instrument object, a dictionary or a filename string')
+ # If we have just loaded a YAML file or constructed from a dictionary, need to convert to correct class
+ for elem_nm, classref in zip(['sample', 'chopper_system', 'moderator', 'detector'],
+ [Sample, ChopperSystem, Moderator, Detector]):
+ try:
+ element = getattr(self, elem_nm)
+ if isinstance(element, dict):
+ setattr(self, elem_nm, classref(element))
+ setattr(self, 'has_' + elem_nm, True)
+ except AttributeError:
+ setattr(self, 'has_' + elem_nm, False)
+ if not self.has_chopper_system or not self.has_moderator:
+ raise AttributeError('No chopper system or moderator found in input.')
+ for method in self.__child_methods:
+ setattr(self, method, getattr(self.chopper_system, method))
+ for prop in self.__child_properties:
+ getter = lambda obj, prop=prop: ChopperSystem.__dict__[prop].__get__(self.chopper_system, ChopperSystem)
+ setter = lambda obj, val, prop=prop: ChopperSystem.__dict__[prop].__set__(self.chopper_system, val)
+ setattr(type(self), prop, property(getter, setter))
+ # Now reset default chopper/variant and frequency
+ if chopper or freq:
+ self.setChopper(chopper if chopper else self.getChopper(), freq if freq else self.frequency)
+
+ def setInstrument(self, instrument):
+ self.__dict__.clear()
+ self.__init__(instrument)
+
+ def getFlux(self, Ei_in=None, frequency=None):
+ """ Returns the monochromatic flux estimate in n/cm^2/s """
+ Ei = self.ei if Ei_in is None else Ei_in
+ return self.moderator.getFlux(Ei) * self.chopper_system.getTransmission(Ei, frequency)
+
+ def getMultiRepFlux(self, Ei_in=None, frequency=None):
+ Ei = self.ei if Ei_in is None else Ei_in
+ if frequency:
+ oldfreq = self.frequency
+ self.frequency = frequency
+ return self.getFlux(self.getAllowedEi(Ei), frequency)
+ if frequency:
+ self.frequency = oldfreq
+
+ def getResFlux(self, Etrans=None, Ei_in=None, frequency=None):
+ """ Returns the resolution and flux as a tuple. """
+ return self.getResolution(Etrans, Ei_in, frequency), self.getFlux(Ei_in, frequency)
+
+ def getWidths(self, Ei_in=None, frequency=None):
+ """ Returns the time FWHM of different components for one rep (Ei) in microseconds """
+ return {k: np.sqrt(v)*1e6 for k, v in list(self.getVanVar(Ei_in, frequency)[1].items())}
+
+ def getMultiWidths(self, Ei_in=None, frequency=None):
+ """ Returns the time FWHM of different components for each possible rep (Ei) in seconds"""
+ Ei = self.ei if Ei_in is None else Ei_in
+ Eis = self.getAllowedEi(Ei)
+ outdic = {'Eis': Eis}
+ widths = [self.getWidths(ei, frequency) for ei in Eis]
+ for ky in widths[0]:
+ outdic[ky] = np.hstack(np.array([w[ky] for w in widths]))
+ return outdic
+
+ def getResolution(self, Etrans=None, Ei_in=None, frequency=None):
+ """
+ Calculates resolution (energy) widths
+
+ van = getResolution()
+ van = getResolution(etrans)
+ van = getResolution(etrans, ei, omega)
+
+ Inputs:
+ etrans - list of numpy array of energy transfers to calculate for (meV) [default: linspace(0.05Ei, 0.95Ei, 19)]
+ ei - incident energy in meV [default: preset energy]
+ omega - chopper frequncy in Hz [default: preset frequency]
+
+ Output:
+ van - the incoherent (Vanadium) energy FWHM at etrans in meV
+ """
+ Ei = self.ei if Ei_in is None else Ei_in
+ # If not set, sets energy transfers to values to compare exactly to RAE's original implementation.
+ if Etrans is None:
+ Etrans = np.linspace(0.05*Ei, 0.95*Ei+0.05*0.05*Ei, 19, endpoint=True)
+ v_van, _ = self.getVanVar(Ei, frequency, Etrans)
+ x2 = self.chopper_system.sam_det
+ Ef = Ei - np.array(Etrans)
+ van = (2 * E2V * np.sqrt(Ef**3 * v_van)) / x2
+ return van
+
+ def getMultiRepResolution(self, Etrans=None, Ei_in=None, frequency=None):
+ """ Returns a list of FWHM in meV for all allowed Ei's in multirep mode (in same order as getAllowedEi)
+ The input energy transfer is interpreted as fractions of Ei. e.g. linspace(0,0.9,100) """
+ if Etrans is None:
+ Etrans = np.linspace(0.05, 0.95, 19, endpoint=True)
+ Ei = self.ei if Ei_in is None else Ei_in
+ return [self.getResolution(Etrans * ei, ei, frequency) for ei in self.getAllowedEi(Ei)]
+
+ def getVanVar(self, Ei_in=None, frequency=None, Etrans=0):
+ """ Calculates the time squared FWHM in s^2 at the sample (Vanadium widths) for different components """
+ Ei = self.ei if Ei_in is None else Ei_in
+ Etrans = np.array(Etrans if np.shape(Etrans) else [Etrans])
+ if frequency and frequency != self.frequency:
+ oldfreq = self.frequency
+ self.frequency = frequency
+ tsqmod = self.moderator.getWidthSquared(Ei)
+ tsqchp = self.chopper_system.getWidthSquared(Ei)
+ tsqjit = self.tjit**2
+ # Gets distances: x0=mod-final chopper, xa=aperture-final, x1=final-sample, x2=sample-det, xm=mod-first chopper
+ x0, xa, x1, x2, xm = self.chopper_system.getDistances()
+ # For Disk chopper spectrometers, the opening times of the first chopper can be the effective moderator time
+ if tsqchp[1] is not None:
+ frac_dist = 1 - (xm / x0)
+ tsmeff = tsqmod * frac_dist**2 # Effective moderator time at first chopper
+ x0 -= xm # Propagate from first chopper, not from moderator (after rescaling tmod)
+ tsqmod = tsmeff if (tsqchp[1] > tsmeff) else tsqchp[1]
+ tsqchp = tsqchp[0]
+ # Propagate the time widths to the detector position
+ omega = self.frequency[0] * 2 * np.pi
+ vi = E2V * np.sqrt(Ei)
+ vf = E2V * np.sqrt(Ei - Etrans)
+ vratio = (vi / vf)**3
+ tanthm = np.tan(self.moderator.theta_m * np.pi / 180.)
+ g1, g2 = tuple(1. - ((omega * tanthm / vi) * np.array([xa + x1, x0 - xa])))
+ f1, f2 = tuple(1. + (x1 / x0) * np.array([g1, g2]))
+ g1, g2, f1, f2 = tuple(np.array([g1, g2, f1, f2]) / (omega * (xa + x1)))
+ modfac = (x1 + vratio * x2) / x0
+ chpfac = 1. + modfac
+ apefac = f1 + ((vratio * x2 / x0) * g1)
+ tsqmod *= modfac**2
+ tsqchp *= chpfac**2
+ tsqjit *= chpfac**2
+ tsqape = apefac**2 * (self.aperture_width**2 / 12.) * SIGMA2FWHMSQ
+ vsqvan = tsqmod + tsqchp + tsqjit + tsqape
+ outdic = {'moderator': tsqmod, 'chopper': tsqchp, 'jitter': tsqjit, 'aperture': tsqape}
+ if self.has_detector:
+ phi = self.detector.phi_deg * np.pi / 180.
+ tsqdet = (1. / vf)**2 * np.array([self.detector.getWidthSquared(Ei, en) for en in Etrans])
+ vsqvan += tsqdet
+ outdic['detector'] = tsqdet
+ else:
+ phi = 0.
+ if self.has_sample:
+ gam = self.sample.gamma_deg * np.pi / 180.
+ bb = (-np.sin(gam) / vi) + (np.sin(gam - phi) / vf) - (f2 * np.cos(gam))
+ samfac = bb - ((vratio * x2 / x0) * g2 * np.cos(gam))
+ tsqsam = samfac**2 * self.sample.getWidthSquared()
+ vsqvan += tsqsam
+ outdic['sample'] = tsqsam
+ if frequency:
+ self.frequency = oldfreq
+ return vsqvan, outdic
+
+ @property
+ def aperture_width(self):
+ if hasattr(self.chopper_system, 'aperture_width') and self.chopper_system.aperture_width:
+ return self.chopper_system.aperture_width
+ return 0.
+
+ @property
+ def aperture_height(self):
+ if hasattr(self.chopper_system, 'aperture_height') and self.chopper_system.aperture_height:
+ return self.chopper_system.aperture_height
+ return 0.
+
+ @property
+ def instname(self):
+ return self.name
+
+ @classmethod
+ def calculate(cls, *args, **kwargs):
+ """
+ ! Calculates the resolution and flux directly (without setting up a PyChop2 object)
+ !
+ ! PyChop2.calculate('mari', 's', 250., 55.) # Instname, Chopper Type, Freq, Ei in order
+ ! PyChop2.calculate('let', 180, 2.2) # For LET, chopper type is not needed.
+ ! PyChop2.calculate('let', [160., 80.], 1.) # For LET, specify resolution and pulse remover freq
+ ! PyChop2.calculate('let', 'High flux', 80, 2.2) # LET default is medium flux configuration
+ ! PyChop2.calculate(inst='mari', package='s', freq=250., ei=55.) # With keyword arguments
+ ! PyChop2.calculate(inst='let', variant='High resolution', freq=[160., 80.], ei=2.2)
+ !
+ ! For LET, the allowed variant names are:
+ ! 'High resolution'
+ ! 'High flux'
+ ! 'Intermediate'
+ ! You have to use these strings exactly.
+ !
+ ! By default this function returns the elastic resolution and flux only.
+ ! If you want the inelastic resolution, specify the inelastic energy transfer
+ ! as either the last positional argument, or as a keyword argument, e.g.:
+ !
+ ! PyChop2.calculate('merlin', 'g', 450., 60., range(55))
+ ! PyChop2.calculate('maps', 'a', 450., 600., etrans=np.linspace(0,550,55))
+ !
+ ! The results are returned as tuple: (resolution, flux)
+ """
+ argdict = argparser(args, kwargs, ['inst', 'package', 'frequency', 'ei', 'etrans', 'variant'])
+ if argdict['inst'] is None:
+ raise RuntimeError('You must specify the instrument name')
+ obj = cls(argdict['inst'])
+ obj.setChopper(argdict['package'], argdict['frequency'])
+ obj.ei = argdict['ei']
+ if argdict['variant']:
+ obj.variant = argdict['variant']
+ return obj.getResolution(argdict['etrans'] if argdict['etrans'] else 0.), obj.getFlux()
+
+ def __repr__(self):
+ return self.name if self.name else 'Undefined instrument'
diff --git a/scripts/PyChop/MulpyRep.py b/scripts/PyChop/MulpyRep.py
index efda669ad6b..5cdedea031f 100644
--- a/scripts/PyChop/MulpyRep.py
+++ b/scripts/PyChop/MulpyRep.py
@@ -184,19 +184,24 @@ def calcChopTimes(efocus, freq, instrumentpars, chop2Phase=5):
chop_times = [] # empty list to hold each chopper opening period
+ # figures out phase information
+ if not (hasattr(ph_ind_v, '__len__') and len(ph_ind_v) == len(dist)):
+ ph_ind = [i==(ph_ind_v[-1] if hasattr(ph_ind_v, '__len__') else ph_ind_v) for i in range(len(dist))]
+ chop2Phase = phase = chop2Phase if hasattr(chop2Phase, '__len__') else [chop2Phase]
+ if len(chop2Phase) != len(dist):
+ if len(chop2Phase) == len(ph_ind_v):
+ phase = [False] * len(dist)
+ for i in range(len(ph_ind_v)):
+ phase[ph_ind_v[i]] = chop2Phase[i]
+ else:
+ phase = [chop2Phase[-1] if ph_ind[i] else False for i in range(len(dist))]
+
# first we optimise on the main Ei
for i in range(len(dist)):
- if hasattr(ph_ind_v, '__len__') and len(ph_ind_v) == len(dist):
- ph_ind = ph_ind_v[i]
- else:
- ph_ind = ph_ind_v if not hasattr(ph_ind_v, '__len__') else ph_ind_v[-1]
# loop over each chopper
- # chopper2 is a special case
- if isinstance(ph_ind, string_types) and i == int(ph_ind) and chop2Phase is not None:
- islt = int(chop2Phase)
- else:
- islt = 0
- if i == ph_ind and chop2Phase is not None:
+ # checks whether this chopper should have an independently set phase / delay
+ islt = int(phase[i]) if (ph_ind[i] and isinstance(phase[i], string_types)) else 0
+ if ph_ind[i] and not isinstance(phase[i], string_types):
# effective chopper velocity (if 2 disks effective velocity is double)
chopVel = 2*np.pi*radius[i] * numDisk[i] * freq[i]
# full opening time
diff --git a/scripts/PyChop/PyChop2.py b/scripts/PyChop/PyChop2.py
index 94b72016269..29bde1f06a7 100644
--- a/scripts/PyChop/PyChop2.py
+++ b/scripts/PyChop/PyChop2.py
@@ -8,6 +8,7 @@ direct geometry time-of-flight inelastic neutron spectrometers.
from __future__ import (absolute_import, division, print_function)
from .ISISFermi import ISISFermi
from .ISISDisk import ISISDisk
+import warnings
class PyChop2:
@@ -35,6 +36,8 @@ class PyChop2:
'MARI': ISISDisk}
def __init__(self, instname, *args):
+ warnings.warn("The PyChop2 class is deprecated and will be removed in the next Mantid version. "
+ "Please use the Instrument class or the official PyChop CLI interface.", DeprecationWarning)
instname = instname.upper()
if instname not in self.__Classes.keys():
raise ValueError('Instrument %s not recognised' % (instname))
diff --git a/scripts/PyChop/PyChopGui.py b/scripts/PyChop/PyChopGui.py
index a2a1bbe1d01..3c9da5650a7 100755
--- a/scripts/PyChop/PyChopGui.py
+++ b/scripts/PyChop/PyChopGui.py
@@ -12,7 +12,7 @@ from __future__ import (absolute_import, division, print_function)
import sys
import re
import numpy as np
-from .PyChop2 import PyChop2
+from .Instruments import Instrument
from PyQt4 import QtGui, QtCore
from matplotlib.backends.backend_qt4agg import FigureCanvasQTAgg as FigureCanvas
from matplotlib.backends.backend_qt4agg import NavigationToolbar2QT as NavigationToolbar
@@ -26,7 +26,7 @@ class PyChopGui(QtGui.QMainWindow):
at spallation sources by calculating the resolution and flux at a given neutron energies.
"""
- instruments = ['MAPS', 'MARI', 'MERLIN', 'LET']
+ instruments = {'MAPS': 'maps.yaml', 'MARI': 'mari.yaml', 'MERLIN': 'merlin.yaml', 'LET': 'let.yaml'}
choppers = {
'MAPS':['A', 'B', 'S'],
'MARI':['A', 'B', 'R', 'G', 'S'],
@@ -95,7 +95,7 @@ class PyChopGui(QtGui.QMainWindow):
else:
self.widgets['MultiRepCheck'].setEnabled(False)
self.widgets['MultiRepCheck'].setChecked(False)
- self.engine.setInstrument(str(instname))
+ self.engine.setInstrument(self.folder + self.instruments[str(instname)])
self.engine.setChopper(str(self.widgets['ChopperCombo']['Combo'].currentText()))
self.engine.setFrequency(float(self.widgets['FrequencyCombo']['Combo'].currentText()))
val = self.flxslder.val * self.maxE[self.engine.instname] / 100
@@ -134,17 +134,17 @@ class PyChopGui(QtGui.QMainWindow):
if 'LET' in str(self.engine.instname):
freqpr = float(self.widgets['PulseRemoverCombo']['Combo'].currentText())
chop2phase = float(self.widgets['Chopper2Phase']['Edit'].text()) % 1e5
- self.engine.setFrequency([freq_in, freqpr], Chopper2Phase=chop2phase)
+ self.engine.setFrequency([freq_in, freqpr], phase=chop2phase)
elif 'MERLIN' in str(self.engine.instname) and 'G' in self.engine.getChopper():
chop2phase = float(self.widgets['Chopper2Phase']['Edit'].text()) % 2e4
- self.engine.setFrequency(freq_in, Chopper2Phase=chop2phase)
+ self.engine.setFrequency(freq_in, phase=chop2phase)
elif 'MAPS' in str(self.engine.instname):
freqpr = float(self.widgets['PulseRemoverCombo']['Combo'].currentText())
- chop2phase = int(self.widgets['Chopper2Phase']['Edit'].text())
- self.engine.setFrequency([freq_in, freqpr], Chopper2Phase=chop2phase)
+ chop2phase = str(self.widgets['Chopper2Phase']['Edit'].text())
+ self.engine.setFrequency([freq_in, freqpr], phase=chop2phase)
elif 'MARI' in str(self.engine.instname):
- chop2phase = int(self.widgets['Chopper2Phase']['Edit'].text())
- self.engine.setFrequency(freq_in, Chopper2Phase=chop2phase)
+ chop2phase = str(self.widgets['Chopper2Phase']['Edit'].text())
+ self.engine.setFrequency(freq_in, phase=chop2phase)
else:
self.engine.setFrequency(freq_in)
@@ -772,7 +772,9 @@ class PyChopGui(QtGui.QMainWindow):
def __init__(self):
super(PyChopGui, self).__init__()
- self.engine = PyChop2('MAPS', 'A', 50)
+ import sys, os
+ self.folder = os.path.dirname(sys.modules[self.__module__].__file__) + '/'
+ self.engine = Instrument(self.folder+self.instruments['MAPS'], 'A', 50)
self.drawLayout()
self.setInstrument('MAPS')
self.xlim = 0
diff --git a/scripts/PyChop/__init__.py b/scripts/PyChop/__init__.py
index 15798b1b95c..e694aedeb47 100644
--- a/scripts/PyChop/__init__.py
+++ b/scripts/PyChop/__init__.py
@@ -19,7 +19,7 @@ PyChop2.showGUI()
from __future__ import (absolute_import, division, print_function)
import warnings
-from .PyChop2 import PyChop2 # noqa: F401
+from .Instruments import Instrument as PyChop2 # noqa: F401
# If the system doesn't have matplotlib, don't import the GUI.
try:
from .PyChopGui import show as showGUI
diff --git a/scripts/PyChop/let.yaml b/scripts/PyChop/let.yaml
new file mode 100644
index 00000000000..5b0bed2544e
--- /dev/null
+++ b/scripts/PyChop/let.yaml
@@ -0,0 +1,138 @@
+name: LET
+# Input file for PyChop2 for the LET spectrometer at ISIS.
+# At a minium, the "chopper_system" and "moderator" entries must be present
+
+chopper_system:
+ name: LET chopper system
+ chop_sam: 1.5 # Distance (x1) from final chopper to sample (m)
+ sam_det: 3.5 # Distance (x2) from sample to detector (m)
+ choppers:
+ - # Each entry must have a dash on an otherwise empty line!
+ name: LET Disk 1 Resolution
+ distance: 7.83 # Distance from moderator to this chopper in metres
+ slot_width: 40 # Slot width in mm
+ guide_width: 40 # Width of guide after chopper in mm
+ nslot: 6 # Number of slots. (Assume all slots equally spaced)
+ radius: 290 # Disk radius
+ isDouble: True # Is this a double disk system?
+ isPhaseIndependent: False # Is this disk to be phased independently?
+ -
+ name: LET Disk 2 Frame Overlap
+ distance: 8.4 # Distance from moderator to this chopper in metres
+ slot_width: 890 # Slot width in mm
+ guide_width: 40 # Width of guide after chopper in mm
+ nslot: 1 # Number of slots. (Assume all slots equally spaced)
+ radius: 545 # Disk radius
+ isDouble: False # Is this a double disk system?
+ isPhaseIndependent: True # Is this disk to be phased independently?
+ defaultPhase: 5 # What is the default phase for this disk (either a time in microseconds
+ # or a slot index for the desired rep to go through).
+ -
+ name: LET Disk 3 Pulse Remover
+ distance: 11.75 # Distance from moderator to this chopper in metres
+ slot_width: 56 # Slot width in mm
+ guide_width: 40 # Width of guide after chopper in mm
+ nslot: 2 # Number of slots. (Assume all slots equally spaced)
+ radius: 290 # Disk radius
+ isDouble: False # Is this a double disk system?
+ isPhaseIndependent: False # Is this disk to be phased independently?
+ -
+ name: LET Disk 4 Contamination Remover
+ distance: 15.66 # Distance from moderator to this chopper in metres
+ slot_width: 56 # Slot width in mm
+ guide_width: 40 # Width of guide after chopper in mm
+ nslot: 6 # Number of slots. (Assume all slots equally spaced)
+ radius: 290 # Disk radius
+ isDouble: False # Is this a double disk system?
+ isPhaseIndependent: False # Is this disk to be phased independently?
+ -
+ name: LET Disk 5 Resolution
+ distance: 23.5 # Distance from moderator to this chopper in metres
+ slot_width: 20 # Slot width in mm
+ guide_width: 20 # Width of guide after chopper in mm
+ nslot: 2 # Number of slots. (Assume all slots equally spaced)
+ radius: 290 # Disk radius
+ isDouble: True # Is this a double disk system?
+ isPhaseIndependent: False # Is this disk to be phased independently?
+ # Now define how the frequencies of the choppers should be related
+ # This is an NxM matrix A where N is the number of choppers and M is the number of indepdent frequencies
+ # Such that A.F will give the N required frequencies for each chopper from the M input frequencies
+ frequency_matrix: # LET is run with two main frequencies:
+ [[0, 0.5], # f1: The frequency of the resolution chopper (Disk 5)
+ [0, 0], # f2: The frequency of the pulse removal chopper (Disk 3)
+ [0, 1], # Disk 4 is usually run at half f1, and disk 1 at half of f2
+ [0.5, 0],
+ [1, 0]]
+ constant_frequencies: # This specifies if a chopper should be run at a fixed constant frequency
+ [0, 10., 0., 0., 0.] # On LET, Disk 2, the frame-overlap chopper, should always run at 10 Hz
+ max_frequencies:
+ [300, 300] # Maximum frequencies (Hz)
+ default_frequencies:
+ [200, 100]
+ overlap_ei_frac: 0.9 # Fraction of energy loss Ei to plot ToF lines in time-distance plots
+ ei_limits: [0, 30] # Limits on ei for multirep calculations (reps outside range ignored)
+ flux_ref_slot: 20 # Reference slot width (mm) for flux transmission calculation (disk choppers only)
+ flux_ref_freq: 150 # Reference final chopper freq (Hz) for flux transmission calc (disk choppers only)
+ # Can define variants which overide one of the above parameters
+ variants: # On LET the final chopper has 3 sets of slots with different widths
+ Intermediate: # This is the default mode (uses the 20mm slot).
+ default: True # If this key is omitted, the variant without any option is used as default
+ High Flux:
+ choppers: # Note that we need to leave the required number of dashes to indicate
+ - # choppers whos parameters have not changed.
+ -
+ -
+ -
+ -
+ slot_width: 31
+ High Resolution:
+ frequency_matrix:
+ [[0.5, 0],
+ [0, 0],
+ [0, 1],
+ [0, 1],
+ [1, 0]]
+ choppers: # Note that we need to leave the required number of dashes to indicate
+ - # choppers whos parameters have not changed.
+ -
+ -
+ -
+ -
+ slot_width: 15
+
+moderator:
+ name: TS2 Hydrogen # A==water, AP==poisoned water, CH4==methane, H2==hydrogen. This is only used for analytical calculations
+ # of the flux distribution for ISIS TS1 moderators. If measured_flux is defined below, name can be anything
+ imod: 3 # Moderator time profile type: 0==chi^2, 1==Ikeda-Carpenter, 2==modified chi^2
+ mod_pars: [-3.143, 49.28, .535] # imod==3 is polynomial. Pars are coeff from highest order to lowest
+ theta: 32.0 # Angle beamline makes with moderator face (degrees)
+ source_rep: 10 # Frequency of source (Hz)
+ emission_time: 3500 # Time width of source pulse in microseconds
+ measured_flux: # Table of measured flux vs wavelength. Wavelength in Angstrom. Flux in n/cm^2/s/uA/Angstrom
+ scale_factor: 63276.8362 # A factor to scale the flux values below by
+ wavelength: [0.5, 0.6, 0.7, 0.8, 0.9, 1. , 1.1, 1.2, 1.3,
+ 1.4, 1.5, 1.6, 1.7, 1.8, 1.9, 2. , 2.1, 2.2,
+ 2.3, 2.4, 2.5, 2.6, 2.7, 2.8, 2.9, 3. , 3.1,
+ 3.2, 3.3, 3.4, 3.5, 3.6, 3.7, 3.8, 3.9, 4. ,
+ 4.1, 4.2, 4.3, 4.4, 4.5, 4.6, 4.7, 4.8, 4.9,
+ 5.0, 5.1, 5.2, 5.3, 5.4, 5.5, 5.6, 5.7, 5.8,
+ 5.9, 6.0, 6.1, 6.2, 6.3, 6.4, 6.5, 6.6, 6.7,
+ 6.8, 6.9, 7.0, 7.1, 7.2, 7.3, 7.4, 7.5, 7.6,
+ 7.7, 7.8, 7.9, 8.0, 8.1, 8.2, 8.3, 8.4, 8.5,
+ 8.6, 8.7, 8.8, 8.9, 9.0, 9.1, 9.2, 9.3, 9.4,
+ 9.5, 9.6, 9.7, 9.8, 9.9, 10.0, 10.1, 10.2, 10.3,
+ 10.4, 10.5, 10.6, 10.7, 10.8, 10.9, 11.0, 11.1, 11.2,
+ 11.3, 11.4, 11.5, 11.6, 11.7, 11.8, 11.9]
+ flux: [0.0889, 0.1003, 0.1125, 0.1213, 0.1274, 0.1358, 0.1455, 0.1562, 0.1702,
+ 0.1902, 0.2149, 0.2496, 0.2938, 0.3537, 0.4315, 0.5244, 0.6415, 0.7856,
+ 0.9341, 1.0551, 1.1437, 1.1955, 1.2004, 1.1903, 1.1662, 1.1428, 1.1176,
+ 1.0875, 1.0641, 1.0562, 1.0242, 0.9876, 0.9586, 0.9415, 0.924, 0.8856,
+ 0.8865, 0.8727, 0.842, 0.8125, 0.7849, 0.7596, 0.7417, 0.7143, 0.6869,
+ 0.6608, 0.6341, 0.6073, 0.581, 0.5548, 0.5304, 0.507, 0.4849, 0.4639,
+ 0.4445, 0.425, 0.407, 0.3902, 0.3737, 0.3579, 0.3427, 0.3274, 0.3129,
+ 0.2989, 0.2854, 0.2724, 0.2601, 0.2483, 0.2371, 0.2267, 0.2167, 0.2072,
+ 0.1984, 0.19, 0.1821, 0.1743, 0.1669, 0.1599, 0.1532, 0.1467, 0.1404,
+ 0.1346, 0.1291, 0.1238, 0.1189, 0.1141, 0.1097, 0.1053, 0.1014, 0.0975,
+ 0.0938, 0.0902, 0.0866, 0.0834, 0.0801, 0.077, 0.0741, 0.0712, 0.0686,
+ 0.066, 0.0637, 0.0614, 0.0593, 0.0571, 0.0551, 0.0532, 0.0512, 0.0494,
+ 0.0477, 0.0461, 0.0445, 0.043, 0.0415, 0.0401, 0.0387]
diff --git a/scripts/PyChop/maps.yaml b/scripts/PyChop/maps.yaml
new file mode 100644
index 00000000000..30b95dbd2ea
--- /dev/null
+++ b/scripts/PyChop/maps.yaml
@@ -0,0 +1,130 @@
+name: MAPS
+# Input file for PyChop2 for the MAPS spectrometer at ISIS.
+
+chopper_system:
+ name: MAPS chopper system
+ chop_sam: 1.9 # Distance (x1) from final chopper to sample (m)
+ sam_det: 6.0 # Distance (x2) from sample to detector (m)
+ aperture_width: 0.094 # Width of aperture at moderator face (m)
+ aperture_height: 0.094 # Height of aperture at moderator face (m)
+ choppers:
+ - # Each entry must have a dash on an otherwise empty line!
+ name: MAPS Disk
+ distance: 8.831 # Distance from moderator to this chopper in metres
+ slot_width: 68 # Slot width in mm
+ guide_width: 50 # Width of guide after chopper in mm
+ nslot: 4 # Number of slots. If slot_ang_pos is specified can ommit this entry
+ # Next line has the angular position of the slots in degrees
+ # Must be monotonically increasing. Can omit if nslot is specifed,
+ # in which case it will be assumed that the slits are equally spaced.
+ slot_ang_pos: [-180., -39.1, 0., 39.1]
+ radius: 375 # Disk radius
+ isDouble: False # Is this a double disk system?
+ isPhaseIndependent: True # Is this disk to be phased independently?
+ defaultPhase: "1" # What is the default phase for this disk (either a time in microseconds
+ # or a slot index [as a string] for the desired rep to go through)
+ -
+ name: MAPS Fermi
+ distance: 10.1 # Distance from moderator to this chopper in metres
+ aperture_distance: 8.27 # Distance from aperture (moderator face) to this chopper (only for Fermi)
+ packages: # A hash of chopper packages
+ A:
+ name: MAPS A (500meV)
+ pslit: 1.087 # Neutron transparent slit width (mm)
+ pslat: 0.534 # Neutron absorbing slat width (mm)
+ radius: 49.0 # Chopper package radius (mm)
+ rho: 1300.0 # Chopper package curvature (mm)
+ tjit: 0.0 # Jitter time (us)
+ fluxcorr: 3.0 # (Empirical/Fudge) factor to scale calculated flux by
+ isPi: False # Should the PI pulse (at 180 deg rotation) be transmitted by this package?
+ B:
+ name: MAPS B (200meV)
+ pslit: 1.812 # Neutron transparent slit width (mm)
+ pslat: 0.534 # Neutron absorbing slat width (mm)
+ radius: 49.0 # Chopper package radius (mm)
+ rho: 920.0 # Chopper package curvature (mm)
+ tjit: 0.0 # Jitter time (us)
+ fluxcorr: 3.0 # (Empirical/Fudge) factor to scale calculated flux by
+ isPi: False # Should the PI pulse (at 180 deg rotation) be transmitted by this package?
+ S:
+ name: MAPS S (Sloppy)
+ pslit: 2.899 # Neutron transparent slit width (mm)
+ pslat: 0.534 # Neutron absorbing slat width (mm)
+ radius: 49.0 # Chopper package radius (mm)
+ rho: 1300.0 # Chopper package curvature (mm)
+ tjit: 0.0 # Jitter time (us)
+ fluxcorr: 2.5 # (Empirical/Fudge) factor to scale calculated flux by
+ isPi: False # Should the PI pulse (at 180 deg rotation) be transmitted by this package?
+ # Now define how the frequencies of the choppers should be related
+ # This is an NxM matrix A where N is the number of choppers and M is the number of indepdent frequencies
+ # Such that A.F will give the N required frequencies for each chopper from the M input frequencies
+ frequency_matrix:
+ [[0, 1], # f1 is the Fermi frequency
+ [1, 0]] # f2 is the Disk frequency
+ max_frequencies:
+ [600, 100] # Maximum frequencies (Hz)
+ default_frequencies:
+ [400, 50]
+
+sample:
+ name: MAPS Sample Can
+ isam: 0 # Sample type: 0==flat plate, 1==ellipse, 2==annulus, 3==sphere, 4==solid cylinder
+ sx: 2.0 # Thickness (mm)
+ sy: 48.0 # Width (mm)
+ sz: 48.0 # Height (mm)
+ gamma: 0.0 # Angle of x-axis to ki (degrees)
+
+detector:
+ name: He3 PSD tubes
+ idet: 2 # Detector type: 1==He tube binned together, 2==He tube
+ dd: 0.025 # Detector depth (diameter for tube) in metres
+ tbin: 0.0 # Detector time bins (microseconds)
+ phi: 0.0 # Detector scattering angle (degrees)
+
+moderator:
+ name: AP # A==water, AP==poisoned water, CH4==methane, H2==hydrogen. This is only used for analytical calculations
+ # of the flux distribution for ISIS TS1 moderators. If measured_flux is defined below, name can be anything
+ imod: 2 # Moderator time profile type: 0==chi^2, 1==Ikeda-Carpenter, 2==modified chi^2
+ mod_pars: [38.6, 0.5226] # Parameters for time profile
+ mod_scale_fn: soft_hat # Function to modify the parameters depending on energy (ommit or leave blank for no modification)
+ mod_scale_par: [1.0, 0.0, 0.0, 150.0, 0.01, 70.0]
+ theta: 32.0 # Angle beamline makes with moderator face (degrees)
+ source_rep: 50 # Frequency of source (Hz)
+ measured_flux: # Table of measured flux vs wavelength. Wavelength in Angstrom. Flux in n/cm^2/s/uA/Angstrom
+ wavelength: [0.0181, 0.0511, 0.0841, 0.1170, 0.1500, 0.1829, 0.2159, 0.2489, 0.2818, 0.3148, 0.3477, 0.3807, 0.4137, 0.4466, 0.4796,
+ 0.5126, 0.5455, 0.5785, 0.6114, 0.6444, 0.6774, 0.7103, 0.7433, 0.7762, 0.8092, 0.8422, 0.8751, 0.9081, 0.9411, 0.9740,
+ 1.0070, 1.0399, 1.0729, 1.1059, 1.1388, 1.1718, 1.2047, 1.2377, 1.2707, 1.3036, 1.3366, 1.3696, 1.4025, 1.4355, 1.4684,
+ 1.5014, 1.5344, 1.5673, 1.6003, 1.6332, 1.6662, 1.6992, 1.7321, 1.7651, 1.7980, 1.8310, 1.8640, 1.8969, 1.9299, 1.9629,
+ 1.9958, 2.0288, 2.0617, 2.0947, 2.1277, 2.1606, 2.1936, 2.2266, 2.2595, 2.2925, 2.3254, 2.3584, 2.3914, 2.4243, 2.4573,
+ 2.4902, 2.5232, 2.5562, 2.5891, 2.6221, 2.6551, 2.6880, 2.7210, 2.7539, 2.7869, 2.8199, 2.8528, 2.8858, 2.9187, 2.9517,
+ 2.9847, 3.0176, 3.0506, 3.0835, 3.1165, 3.1495, 3.1824, 3.2154, 3.2484, 3.2813, 3.3143, 3.3472, 3.3802, 3.4132, 3.4461,
+ 3.4791, 3.5120, 3.5450, 3.5780, 3.6109, 3.6439, 3.6769, 3.7098, 3.7428, 3.7757, 3.8087, 3.8417, 3.8746, 3.9076, 3.9406,
+ 3.9735, 4.0065, 4.0394, 4.0724, 4.1054, 4.1383, 4.1713, 4.2042, 4.2372, 4.2702, 4.3031, 4.3361, 4.3690, 4.4020, 4.4350,
+ 4.4679, 4.5009, 4.5339, 4.5668, 4.5998, 4.6327, 4.6657, 4.6987, 4.7316, 4.7646, 4.7976, 4.8305, 4.8635, 4.8964, 4.9294,
+ 4.9624, 4.9953, 5.0283, 5.0612, 5.0942, 5.1272, 5.1601, 5.1931, 5.2260, 5.2590, 5.2920, 5.3249, 5.3579, 5.3909, 5.4238,
+ 5.4568, 5.4897, 5.5227, 5.5557, 5.5886, 5.6216, 5.6546, 5.6875, 5.7205, 5.7534, 5.7864, 5.8194, 5.8523, 5.8853, 5.9182,
+ 5.9512, 5.9842, 6.0171, 6.0501, 6.0831, 6.1160, 6.1490, 6.1819, 6.2149, 6.2479, 6.2808, 6.3138, 6.3467, 6.3797, 6.4127,
+ 6.4456, 6.4786, 6.5115, 6.5445, 6.5750]
+ flux: [2.60775e+08, 8.08616e+07, 4.59455e+07, 3.16183e+07, 2.37830e+07, 1.91458e+07, 1.58204e+07, 1.36088e+07, 1.19225e+07,
+ 1.06105e+07, 9.55533e+06, 8.77411e+06, 8.21901e+06, 7.67656e+06, 7.15055e+06, 6.70757e+06, 6.39276e+06, 6.22826e+06,
+ 6.19798e+06, 6.27731e+06, 6.47994e+06, 6.91604e+06, 7.50573e+06, 8.09035e+06, 8.76875e+06, 9.49975e+06, 1.01658e+07,
+ 1.07548e+07, 1.13597e+07, 1.19941e+07, 1.25374e+07, 1.28821e+07, 1.31248e+07, 1.32727e+07, 1.32305e+07, 1.30834e+07,
+ 1.27595e+07, 1.24351e+07, 1.20115e+07, 1.15789e+07, 1.11218e+07, 1.07191e+07, 1.03272e+07, 9.93856e+06, 9.59376e+06,
+ 9.19836e+06, 8.81571e+06, 8.50457e+06, 8.29274e+06, 8.05058e+06, 7.88437e+06, 7.63907e+06, 7.32047e+06, 6.99681e+06,
+ 6.68968e+06, 6.45270e+06, 6.24161e+06, 6.01323e+06, 5.74713e+06, 5.48816e+06, 5.27153e+06, 5.01780e+06, 4.76127e+06,
+ 4.48172e+06, 4.21345e+06, 3.97093e+06, 3.74819e+06, 3.53683e+06, 3.32935e+06, 3.12404e+06, 2.92801e+06, 2.74479e+06,
+ 2.61634e+06, 2.48606e+06, 2.38826e+06, 2.29410e+06, 2.17636e+06, 2.07461e+06, 1.97063e+06, 1.87220e+06, 1.77780e+06,
+ 1.70202e+06, 1.62584e+06, 1.55763e+06, 1.48989e+06, 1.42924e+06, 1.37959e+06, 1.34056e+06, 1.31926e+06, 1.28573e+06,
+ 1.25559e+06, 1.22426e+06, 1.18988e+06, 1.15714e+06, 1.13032e+06, 1.09423e+06, 1.06161e+06, 1.02650e+06, 9.95519e+05,
+ 9.56437e+05, 9.24815e+05, 8.90446e+05, 8.56656e+05, 8.28196e+05, 8.01094e+05, 7.79358e+05, 7.56306e+05, 7.35949e+05,
+ 7.24375e+05, 7.02174e+05, 6.86458e+05, 6.65894e+05, 6.43176e+05, 6.24539e+05, 6.01304e+05, 5.82505e+05, 5.61653e+05,
+ 5.41996e+05, 5.25903e+05, 5.10613e+05, 4.96677e+05, 4.82118e+05, 4.64661e+05, 4.65809e+05, 4.68617e+05, 4.56137e+05,
+ 4.42141e+05, 4.27460e+05, 4.10041e+05, 3.98628e+05, 3.84161e+05, 3.71166e+05, 3.57501e+05, 3.45980e+05, 3.35925e+05,
+ 3.23733e+05, 3.13815e+05, 3.03413e+05, 2.91757e+05, 2.82348e+05, 2.72917e+05, 2.57271e+05, 2.41863e+05, 2.58619e+05,
+ 2.53316e+05, 2.43464e+05, 2.35779e+05, 2.29787e+05, 2.22481e+05, 2.14144e+05, 2.08181e+05, 2.01866e+05, 1.95864e+05,
+ 1.89808e+05, 1.84740e+05, 1.78016e+05, 1.73397e+05, 1.68777e+05, 1.65549e+05, 1.61071e+05, 1.56594e+05, 1.52364e+05,
+ 1.49546e+05, 1.46162e+05, 1.43155e+05, 1.40167e+05, 1.37583e+05, 1.35294e+05, 1.32192e+05, 1.30639e+05, 1.27633e+05,
+ 1.25179e+05, 1.23187e+05, 1.21203e+05, 1.18074e+05, 1.15095e+05, 1.12187e+05, 1.10561e+05, 1.08411e+05, 1.05109e+05,
+ 1.03695e+05, 1.01165e+05, 9.87797e+04, 9.77841e+04, 9.40768e+04, 9.27353e+04, 9.14937e+04, 8.88289e+04, 8.74353e+04,
+ 8.53251e+04, 8.30339e+04, 8.22249e+04, 7.94099e+04, 7.79037e+04, 7.62865e+04, 7.47047e+04, 7.38535e+04, 7.17228e+04,
+ 6.98927e+04, 6.89509e+04]
--
GitLab