Skip to content

Commit

Permalink
Merge pull request #81 from ULB-Metronu/issue-78
Browse files Browse the repository at this point in the history
  • Loading branch information
Robin Tesse committed Dec 31, 2023
2 parents 785bb2a + 3095b9c commit 46ab78c
Show file tree
Hide file tree
Showing 8 changed files with 692 additions and 31 deletions.
8 changes: 5 additions & 3 deletions docs/changelog.rst
View file
Original file line number Diff line number Diff line change
Expand Up @@ -4,19 +4,21 @@ Changelog

* 2023-2

* Add an example with energy degradation.
* Fix filename to compute the coefficients with BDSIM.
* Improve the documentation for the apertures
* Add the paper published in Software-X
* Update librairies

.. todolist::

* 2023-1

* Add documentation of the PTW module
* Add documentation for the optimisation of a beamline
* Minor fix in observers
* Add a plotting method for the SymmetryObserver

.. todolist::


* 2022-1

* Write the documentation
Expand Down
7 changes: 4 additions & 3 deletions docs/modules/manzoni.rst
View file
Original file line number Diff line number Diff line change
Expand Up @@ -32,13 +32,14 @@ The module is structured as depicted in figure below and described in the next s
manzoni/structure


Example
=======
Examples
========

.. toctree::
:maxdepth: 1

manzoni/manzoni_example
manzoni/example_tracking
manzoni/example_energy_degradation

Validation
==========
Expand Down
206 changes: 206 additions & 0 deletions docs/modules/manzoni/example_energy_degradation.rst
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,206 @@
Particle tracking with energy degradation
#########################################

We replace the collimator from the previous example by an 10cm length energy degrader in Beryllium.
The parameters `WITH_LOSSES` can be set to `True`` or `False`` depending the application. If it is set to `True`,
the particles are randomly lost. The number of lost particles depends on the material and the kinetic energy and
the coefficient are computed as described in the Fermi module.

.. jupyter-input::

c1 = georges.Element.Degrader(NAME='DEG',
MATERIAL=georges.fermi.materials.Beryllium,
L = 10*_ureg.cm,
WITH_LOSSES=True
)

Let's first import the necessary packages

.. jupyter-execute::
:hide-output:

import matplotlib.pyplot as plt

import georges
from georges import ureg as _ureg
from georges import vis
from georges.manzoni import Input, observers
from georges.manzoni.beam import MadXBeam

Let's define the line using a :code:`PlacementSequence`
-------------------------------------------------------

.. jupyter-execute::
:hide-output:

d1 = georges.Element.Drift(
NAME="D1",
L=0.3 * _ureg.m,
APERTYPE="RECTANGULAR",
APERTURE=[5 * _ureg.cm, 3 * _ureg.cm],
)

qf = georges.Element.Quadrupole(
NAME="Q1",
L=0.3 * _ureg.m,
K1=2 * _ureg.m**-2,
APERTYPE="RECTANGULAR",
APERTURE=[5 * _ureg.cm, 3 * _ureg.cm],
)

d2 = georges.Element.Drift(
NAME="D2",
L=0.3 * _ureg.m,
APERTYPE="CIRCULAR",
APERTURE=[5 * _ureg.cm, 5 * _ureg.cm],
)

b1 = georges.Element.SBend(
NAME="B1",
L=1 * _ureg.m,
ANGLE=30 * _ureg.degrees,
K1=0 * _ureg.m**-2,
APERTYPE="CIRCULAR",
APERTURE=[5 * _ureg.cm, 5 * _ureg.cm],
)

d3 = georges.Element.Drift(
NAME="D3",
L=0.3 * _ureg.m,
APERTYPE="CIRCULAR",
APERTURE=[5 * _ureg.cm, 5 * _ureg.cm],
)

qd = georges.Element.Quadrupole(
NAME="Q2",
L=0.3 * _ureg.m,
K1=-2 * _ureg.m**-2,
APERTYPE="RECTANGULAR",
APERTURE=[5 * _ureg.cm, 3 * _ureg.cm],
)

d4 = georges.Element.Drift(
NAME="D4",
L=0.3 * _ureg.m,
APERTYPE="CIRCULAR",
APERTURE=[5 * _ureg.cm, 5 * _ureg.cm],
)

c1 = georges.Element.Degrader(NAME='DEG',
MATERIAL=georges.fermi.materials.Beryllium,
L = 10*_ureg.cm,
WITH_LOSSES=True
)

d5 = georges.Element.Drift(
NAME="D5",
L=0.3 * _ureg.m,
APERTYPE="CIRCULAR",
APERTURE=[5 * _ureg.cm, 5 * _ureg.cm],
)

b2 = georges.Element.SBend(
NAME="B2",
L=1 * _ureg.m,
ANGLE=-30 * _ureg.degrees,
K1=0 * _ureg.m**-2,
APERTYPE="RECTANGULAR",
APERTURE=[5 * _ureg.cm, 3 * _ureg.cm],
)

d6 = georges.Element.Drift(
NAME="D6",
L=0.3 * _ureg.m,
APERTYPE="CIRCULAR",
APERTURE=[5 * _ureg.cm, 5 * _ureg.cm],
)

sequence = georges.PlacementSequence(name="Sequence")

sequence.place(d1, at_entry=0 * _ureg.m)
sequence.place_after_last(qf)
sequence.place_after_last(d2)
sequence.place_after_last(b1)
sequence.place_after_last(d3)
sequence.place_after_last(c1)
sequence.place_after_last(d4)
sequence.place_after_last(qd)
sequence.place_after_last(d5)
sequence.place_after_last(b2)
sequence.place_after_last(d6)

We use a Gaussian beam with an energy of 230 MeV
------------------------------------------------

.. jupyter-execute::
:hide-output:

kin = georges.Kinematics(230 * _ureg.MeV, particle=georges.particles.Proton, kinetic=True)
sequence.metadata.kinematics = kin

beam = MadXBeam(
kinematics=kin,
distribution=georges.Distribution.from_5d_multigaussian_distribution(
n=10000, xrms=0.1 * _ureg.cm, yrms=0.7 * _ureg.cm, pxrms=0.01, pyrms=0.01
).distribution.values,
)

We can now track in our line with :code:`Manzoni`
-------------------------------------------------

.. jupyter-execute::
:hide-output:

mi = Input.from_sequence(sequence=sequence)

We need to adjust the energy of the line in presence of an energy degrader.

.. jupyter-execute::
:hide-output:

mi.adjust_energy(input_energy=kin.ekin)

Let's freeze the sequence and run Manzoni

.. jupyter-execute::
:hide-output:

mi.freeze()
beam_observer_std = mi.track(beam=beam, observers=observers.SigmaObserver())
beam_observer_beam = mi.track(beam=beam, observers=observers.BeamObserver(with_input_beams=True))
beam_observer_losses = mi.track(beam=beam, observers=observers.LossesObserver())

Plot results
------------

.. tabs::

.. tab:: Standard Deviation

.. jupyter-execute::

fig = plt.figure(figsize=(10,4))
ax = fig.add_subplot(111)
manzoni_plot = vis.ManzoniMatplotlibArtist(ax=ax)
manzoni_plot.plot_cartouche(sequence.df)
manzoni_plot.tracking(beam_observer_std, plane="both")

.. tab:: Losses

.. jupyter-execute::

fig = plt.figure(figsize=(10,4))
ax = fig.add_subplot(111)
manzoni_plot = vis.ManzoniMatplotlibArtist(ax=ax)
manzoni_plot.plot_cartouche(sequence.df)
manzoni_plot.losses(beam_observer_losses, log_scale=False)

.. tab:: Phase-space

.. jupyter-execute::

fig = plt.figure(figsize=(10,4))
ax = fig.add_subplot(111)
manzoni_plot = vis.ManzoniMatplotlibArtist(ax=ax)
manzoni_plot.plot_cartouche(sequence.df)
manzoni_plot.phase_space(beam_observer_beam, element="D5")
4 changes: 2 additions & 2 deletions docs/modules/manzoni/manzoni_example.rst → docs/modules/manzoni/example_tracking.rst
View file
Original file line number Diff line number Diff line change
@@ -1,5 +1,5 @@
Basic example
#############
Particle tracking with collimator
#################################

In this example, we define beamline with `Quadrupoles`, `SBEND` and `Collimators` and we use different observers to analyse the beam and its properties.

Expand Down
95 changes: 93 additions & 2 deletions georges/manzoni/elements/collimators.py
View file
Original file line number Diff line number Diff line change
@@ -1,35 +1,126 @@
"""
TODO
"""
from ... import ureg as _ureg
from .magnets import Drift as _Drift


class Collimator(_Drift):
"""
Define a Collimator.
Attributes:
PARAMETERS (dict): Dictionary containing the parameters of the Collimator with their default values.
Examples:
>>> c1 = Collimator('C1', L=10*_ureg.cm, APERTYPE='ELIPTICAL', APERTURE=[2*_ureg.cm, 3*_ureg.cm])
>>> c1 #doctest: +NORMALIZE_WHITESPACE
Collimator: {'NAME': 'C1',
'AT_ENTRY': <Quantity(0, 'meter')>,
'AT_CENTER': <Quantity(0, 'meter')>,
'AT_EXIT': <Quantity(0, 'meter')>,
'L': <Quantity(10, 'centimeter')>,
'APERTYPE': 'ELIPTICAL',
'APERTURE': [<Quantity(2, 'centimeter')>, <Quantity(3, 'centimeter')>]}
"""

PARAMETERS = {
"APERTYPE": (None, "Aperture type (CIRCULAR, ELIPTIC or RECTANGULAR)"),
"APERTYPE": (None, "Aperture type (CIRCULAR, ELLIPTICAL, RECTANGULAR or PHASE_SPACE)"),
"APERTURE": ([None, None, None, None], ""),
}


class CircularCollimator(Collimator):
"""
Define a CircularCollimator.
Attributes:
PARAMETERS (dict): Dictionary containing the parameters of the CircularCollimator with their default values.
Examples:
>>> c2 = CircularCollimator('C2', L=10*_ureg.cm, APERTURE=[2*_ureg.cm])
>>> c2 #doctest: +NORMALIZE_WHITESPACE
CircularCollimator: {'NAME': 'C2',
'AT_ENTRY': <Quantity(0, 'meter')>,
'AT_CENTER': <Quantity(0, 'meter')>,
'AT_EXIT': <Quantity(0, 'meter')>,
'L': <Quantity(10, 'centimeter')>,
'APERTYPE': 'CIRCULAR',
'APERTURE': [<Quantity(2, 'centimeter')>]}
"""

PARAMETERS = {
"APERTYPE": ("CIRCULAR", "Aperture type"),
}


class EllipticalCollimator(Collimator):
"""
Define an EllipticalCollimator.
Attributes:
PARAMETERS (dict): Dictionary containing the parameters of the EllipticalCollimator with their default values.
Examples:
>>> c2 = EllipticalCollimator('C2', L=10*_ureg.cm, APERTURE=[2*_ureg.cm, 3*_ureg.cm])
>>> c2 #doctest: +NORMALIZE_WHITESPACE
EllipticalCollimator: {'NAME': 'C2',
'AT_ENTRY': <Quantity(0, 'meter')>,
'AT_CENTER': <Quantity(0, 'meter')>,
'AT_EXIT': <Quantity(0, 'meter')>,
'L': <Quantity(10, 'centimeter')>,
'APERTYPE': 'ELLIPTICAL',
'APERTURE': [<Quantity(2, 'centimeter')>, <Quantity(3, 'centimeter')>]}
"""

PARAMETERS = {
"APERTYPE": ("ELLIPTICAL", "Aperture type"),
}


class RectangularCollimator(Collimator):
"""
Define a RectangularCollimator.
Attributes:
PARAMETERS (dict): Dictionary containing the parameters of the RectangularCollimator with their default values.
Examples:
>>> c3 = RectangularCollimator('C3', L=10*_ureg.cm, APERTURE=[2*_ureg.cm, 3*_ureg.cm])
>>> c3 #doctest: +NORMALIZE_WHITESPACE
RectangularCollimator: {'NAME': 'C3',
'AT_ENTRY': <Quantity(0, 'meter')>,
'AT_CENTER': <Quantity(0, 'meter')>,
'AT_EXIT': <Quantity(0, 'meter')>,
'L': <Quantity(10, 'centimeter')>,
'APERTYPE': 'RECTANGULAR',
'APERTURE': [<Quantity(2, 'centimeter')>, <Quantity(3, 'centimeter')>]}
"""

PARAMETERS = {
"APERTYPE": ("RECTANGULAR", "Aperture type"),
}


class Dump(RectangularCollimator):
"""
Define a Dump.
Attributes:
PARAMETERS (dict): Dictionary containing the parameters of the Dump with their default values.
Examples:
>>> d1 = Dump('D1', L=10*_ureg.cm)
>>> d1 #doctest: +NORMALIZE_WHITESPACE
Dump: {'NAME': 'D1',
'AT_ENTRY': <Quantity(0, 'meter')>,
'AT_CENTER': <Quantity(0, 'meter')>,
'AT_EXIT': <Quantity(0, 'meter')>,
'L': <Quantity(10, 'centimeter')>,
'APERTYPE': 'RECTANGULAR',
'APERTURE': <Quantity(0.0, 'centimeter')>}
"""

PARAMETERS = {
"APERTURE": [0.0, 0.0],
"APERTURE": [0.0 * _ureg.cm, 0.0 * _ureg.cm],
}
Loading

0 comments on commit 46ab78c

Please sign in to comment.