Commit aa6bda72 authored by Zolnierczuk, Piotr's avatar Zolnierczuk, Piotr
Browse files

flipper improvements

parent 6bf98f2c
Loading
Loading
Loading
Loading
+2 −2
Original line number Diff line number Diff line
@@ -64,8 +64,8 @@ J2_MAX = 0.0200 # [T*m] max field integral (20 mTm = 200 Gauss*m) # TODO: verif
#
J_LIMS = (J2_MIN, J1_MAX)

# 4b. Flipper names
NSE_FLIPPER_NAMES = {
# 4b. standard NSE flippers. Position measured from the sample (negative is upstream, positive downstream)
NSE_FLIPPERS = { #  type,        position [m],        thickness [m],   turns_density [turns/m]
    "fpi21"  : dict(type='pi/2', position =-2.824663, thickness=0.010, turns_density=833.33), #w=0.10 h=0.18, i33
    "fpi21rs": dict(type='pi'  , position =-2.283902, thickness=0.010, turns_density=833.33), #w=0.10 h=0.18, i37
    "fpi"    : dict(type='pi'  , position =-0.146000, thickness=0.013, turns_density=833.33), #w=0.18 h=0.18, i32
+32 −14
Original line number Diff line number Diff line
@@ -7,14 +7,8 @@ import numpy as np
import scipy.optimize as opt
from numpy import sqrt, cos

try:
from ..constants import (PI, MU_0, HMN, OMEGA_N)
except ImportError:
    from scipy import constants as _const
    PI      = np.pi
    MU_0    = 4e-7*PI
    HMN     = _const.Planck/_const.neutron_mass
    OMEGA_N = _const.physical_constants['neutron gyromag. ratio'][0]/HMN
from ..config import NSE_FLIPPERS

PISQRT2 = PI/sqrt(2)

@@ -41,6 +35,14 @@ class BaseFlipper:
        "the flipper"
        return f"{self.name:8s} (pos={self.position:g}m d={self.thickness:g}m n={self.turns_density:g})"

    def __eq__(self, value):
        "compare two flippers"
        assert isinstance(value, BaseFlipper)
        # ignore the name
        return (self.position      == value.position      and
                self.thickness     == value.thickness     and
                self.turns_density == value.turns_density    )

    def current(self, arg, rotation=None):
        """calculate flipper current
            arg      - wavelength or magnetic field (if rotation is None)
@@ -50,10 +52,14 @@ class BaseFlipper:
        bfield = arg if rotation is None else rotation/(OMEGA_N*self.thickness*arg)
        return bfield/current_sheet(self.turns_density)/2 # factor of 2 for two sheets per flipper


class FlipperPi(BaseFlipper):
    """Simple pi-flipper class"""
    def curve(self, lam, lam0, b_ext=None):

    def current(self, arg, rotation=PI):
        "flipper current"
        return super().current(arg, rotation=rotation)

    def curve(self, lam, lam0, **kwargs):
        """pi-flipper ramp curve
        simple scaling by neutron wavelength
        also takes into account external fields
@@ -62,7 +68,7 @@ class FlipperPi(BaseFlipper):
        lam0   - reference wavelength (lambda_max at SNS-NSE)
        bx_ext - external field at the flipper (beam direction, horizontal axis x)
        bz_ext - external field at the flipper (vertical direction)"""

        b_ext = kwargs.pop('b_ext', None)
        (bx_ext, _ , bz_ext) = (0, 0, 0) if b_ext is None else b_ext

        sigma = PI/OMEGA_N/self.thickness
@@ -71,9 +77,6 @@ class FlipperPi(BaseFlipper):

        return fl/f0

    def current(self, arg, rotation=PI):
        "flipper current"
        return super().current(arg, rotation=rotation)

class FlipperPi2(BaseFlipper):
    """Simple pi/2-flipper class """
@@ -92,6 +95,21 @@ class FlipperPi2(BaseFlipper):
        omega_z = np.asarray([opt.newton(rxx13,PI/2,args=(x,)) for x in omega_x])
        return omega_z/omega_x

def make_flipper(name, params=None):
    """Make a flipper from a config dictionary"""
    if not params:
        params = NSE_FLIPPERS.get(name, {})
        if not params:
            raise ValueError(f"unknown flipper name '{name}', choose from {list(NSE_FLIPPERS.keys())}")
    # make a copy, so that we can "remove" type
    params = params.copy()
    ftype = params.pop('type', 'unknown')
    if ftype == 'pi':
        return FlipperPi(name=name, **params)
    if ftype  == 'pi/2':
        return FlipperPi2(name=name, **params)
    raise ValueError(f"unknown flipper type '{ftype}' for flipper '{name}'")

def tof_curve(flipper, lmax, dlam, points, sample_to_source, **kwargs):
    """
    flipper - flipper object
+11 −29
Original line number Diff line number Diff line
@@ -6,18 +6,8 @@ import numpy as np
import matplotlib.pyplot as plt

from ..constants import GAUSS
from ..config import L2, Ltot, HMN, ANGSTROM, NSE_FLIPPER_NAMES
from ..echo.flippers import FlipperPi, FlipperPi2, tof_curve

#NSE_FLIPPERS = {
#    "fpi21"  : FlipperPi2('fpi21',   position =-2.824663, thickness=0.010, turns_density=833.33), #w=0.10 h=0.18, i33
#    "fpi21rs": FlipperPi ('fpi21rs', position =-2.283902, thickness=0.010, turns_density=833.33), #w=0.10 h=0.18, i37
#    "fpi"    : FlipperPi ('fpi',     position =-0.146000, thickness=0.013, turns_density=833.33), #w=0.18 h=0.18, i32
#    "fpi22rs": FlipperPi ('fpi22rs', position =+2.279194, thickness=0.010, turns_density=821.40), #w=0.28 h=0.28, i38
#    "fpi22"  : FlipperPi2('fpi22',   position =+2.823903, thickness=0.010, turns_density=821.40), #w=0.28 h=0.28, i34
#}

NSE_FLIPPERS = {}
from ..config import L2, Ltot, HMN, ANGSTROM
from ..echo.flippers import make_flipper, tof_curve

def plot_flipper_curve(axs, name='fpi', **kwargs):
    "test with plots"
@@ -51,30 +41,21 @@ def plot_flipper_curve(axs, name='fpi', **kwargs):
    delam = HMN/(freq*L) # wavelengtgh band
    lmax = lmax*ANGSTROM

    params = NSE_FLIPPER_NAMES.get(name, {})
    if not params:
        raise ValueError(f"unknown flipper name '{name}', choose from {list(NSE_FLIPPER_NAMES.keys())}")
    flipper_type = params.pop('type', 'unknown')
    if flipper_type == 'pi':
        f = FlipperPi(name=name, **params)
    elif flipper_type  == 'pi/2':
        f = FlipperPi2(name=name, **params)
    else:
        raise ValueError(f"unknown flipper type '{params.get('type')}' for flipper '{name}'")

    flipper = make_flipper(name)

    curve =  tof_curve(f, lmax, delam, points=points, sample_to_source=L1,
    curve =  tof_curve(flipper, lmax, delam, points=points, sample_to_source=L1,
                       freq=freq, b_ext=B_ext, rewind=5)
    curve = np.roll(curve, shift)
    scurve = curve*f.current(lmax)*scale
    scurve = curve*flipper.current(lmax)*scale
    p = axs.plot(scurve,'-', lw=1, label='calc %s' % name)
    if compare is not None:
        axs.plot(compare_data,'--', lw=1, label='compare %s' % compare, color=p[0].get_color())

    if out:
        with open(out, 'wt', encoding='ascii') as fd:
            fd.write(f"# flipper {f}\n")
            fd.write(f"# lambda = ({lmax-delam/ANGSTROM:g} - {lmax:g}) A\n")
            fd.write(f"# flipper {flipper}\n")
            fd.write(f"# instrument {pos} (L={L:g}m)\n")
            fd.write(f"# lambda ({(lmax-delam)/ANGSTROM:g} - {lmax/ANGSTROM:g}) A\n")
            for point in scurve:
                fd.write(f"{point:g}\n")

@@ -87,9 +68,10 @@ def plot_flipper_curve(axs, name='fpi', **kwargs):
def action_flipper(_filenames, **kwargs):
    "default action for echo, atari, etc"
    savefig  = kwargs.pop('savefig')
    savefile = kwargs.pop('savefile')

    fig, ax = plt.subplots(1,1,figsize=(8,6))
    title, _ = plot_flipper_curve(ax, **kwargs)
    title, _ = plot_flipper_curve(ax, out=savefile, **kwargs)
    fig.suptitle(title)

    if savefig:
+1 −1
Original line number Diff line number Diff line
@@ -50,7 +50,7 @@ def norm_pix(x1, x2, nx, whole_detector=False):
    "normalize pixel"
    if whole_detector:
        return 0, nx
    return np.mod(x1,nx), np.mod(x2,nx)
    return np.mod(x1,nx+1), np.mod(x2,nx+1)

def get_pix(kwargs):
    """get default pixels, TOF bins"""
+2 −2
Original line number Diff line number Diff line
@@ -3,8 +3,8 @@ PySEN revision module
"""
import sys
__version__  = "2.0"
__release__  = ".7b3"
__date__     = "Oct 10, 2025"
__release__  = ".7b4"
__date__     = "Oct 14, 2025"

def version(full=False):
    "get pysen version number"
Loading