Bug fix in atgeometry #667

lfarv · 2023-10-05T09:42:46Z

Fixes a bug in atgeometry reported by @simoneliuzzo.

The centered flag now works as expected.

However, the computation of the radius is made using the 1st and last points. For a full ring, these points almost coincide, which makes the computation very imprecise. This happens for both Matlab and python.

For better results, the computation of the radius should be completely rewritten by taking all points into account (noting that this "polygonal" radius is not a constant…). This if left for another pull request.

simoneliuzzo · 2023-10-09T07:03:12Z

The centered option seems to work a bit better, but not yet as expected.

lcarver · 2023-10-09T07:26:12Z

atmat/lattice/survey/atgeometry.m

 if centered
-    yy=yy+radius;
+    xx=xx-radius*sin(theta0);


maybe i'm missing something. But why not do

xcenter = (max(xx) + min(xx) )/2; ycenter = (max(yy) + min(yy) )/2; xx=xx-xcenter; yy=yy-ycenter;

@lcarver would work only for full SR.

simoneliuzzo · 2023-10-10T10:49:56Z

Dear @lfarv,

the centered option seems to 'center' on ~2.5*10^5 rather than on ~140.

may be you are already looking at this issue.

best regards
simone

lfarv · 2023-10-10T12:20:46Z

@simoneliuzzo : see my previous answer: the present geometric computation cannot work for full lattices. So

we stay as now: this bug fix gives at least a decent plot, and you can still use offset to put the center where you like,
we make a special case for full lattices with one of the methods given above.

simoneliuzzo · 2023-10-10T12:27:04Z

Dear @lfarv,
I have memory of this feature working correctly. In pyAT for example it is working as expected. Could we apply the same strategy?
regards
Simone

lfarv · 2023-10-10T14:38:07Z

pyAT suffers from exactly the same problem !

As I already wrote, for a full ring, you cannot determine the center by looking at both ends.

The radius being defined as the intersection of the perpendiculars to both end of the lattice, for a full ring it goes to infinity. If because of rounding errors the ring does not close exactly, then you get any crazy value for the radius, it depends on your ring…

If we go for a special case for full rings, it will of course be applied to both Matlab and python.

swhite2401 · 2023-10-11T08:43:37Z

Since we are looking into this it would probably make sense to fix it here.

I could not find the proposed solutions mentioned in a previous message but since the radius is not constant I suppose we would have to compute some kind of average, was this the proposed solution? This could be quite heavy calculation no?

carmignani · 2023-10-11T11:45:10Z

I think what @lcarver proposed is the best: the centre would be the average of max and min xx and yy positions.

lfarv · 2023-10-11T16:33:58Z

My proposal: if the input is a cell of a periodic ring, the centre is non-ambiguous and the present solution is correct, we keep it. If fails in 2 cases: full ring and half-ring, when the perpendiculars to both ends coincide. So we need 2 special cases:

Full ring: theta-theta0 ≈ 2pi. There is no well-defined centre, we take either Lee's solution or take the average of x values and y values. This is good enough.
Half ring: theta-theta0 ≈ pi. An approximate center is the middle of both ends.

I can't help until next week, so if anyone wants to deal with that, he's welcome !

lfarv · 2023-10-23T15:32:31Z

Finally, centring now works in most cases, including full rings and half rings. The fix is applied to both Matlab and python.

It's now ready for merging. @lcarver and @simoneliuzzo, any comment?

swhite2401 · 2023-10-24T10:40:39Z

I looked at the python and this seems fine, I just have one cosmetic comment _EPSIL is defined a the top level while it will be used only locally in the function.
I can imagine that other function would need an epsilon in the future so defining it at the top level may not be the best.... for now this is fine anyway so it is up to you.

lfarv · 2023-10-24T15:48:56Z

@swhite2401: defining the threshold globally makes it visible in case one wants to modify it, and makes it assigned once at import rather that each time the function is called.

To avoid the ambiguity that you pointed out, I renamed it to _GEOMETRY_EPSIL.

Waiting for approval ! @simoneliuzzo ?

simoneliuzzo · 2023-10-25T11:41:32Z

The python version gives me this error:

/usr/local/bin/python3.9 /Users/liuzzo/PycharmProjects/testatgeometry/testgeom.py 
Traceback (most recent call last):
  File "/Users/liuzzo/PycharmProjects/testatgeometry/testgeom.py", line 1, in <module>
    import at
  File "/Users/liuzzo/beamdyn/at/pyat/at/__init__.py", line 3, in <module>
    from ._version import __version__, __version_tuple__
ModuleNotFoundError: No module named 'at._version'

I installed with these commands:
python3.9 -m venv atgeo
source atgeo/bin/activate
pip install --upgrade pip
pip install git+https://github.com/atcollab/at.git@fix_atgeometry
pip install matplotlib

I will try the matlab version now.

simoneliuzzo · 2023-10-25T11:52:13Z

The matlab version gives a circle for centered but vertically displaced by 2.5 10^5 (as in the previous version)

to install I did

git clone https://github.com/atcollab/at.git at_geom_25Oct
git checkout fix_atgeometry

>> restoredefaultpath;
>> addpath(genpath('./at_geom_25Oct/atmat'))
>> addpath(genpath('./at_geom_25Oct/atintegrators'))
>> atmexall;

lfarv · 2023-10-25T15:04:08Z

@simoneliuzzo: for your Matlab test: apparently you are running an old version. The git checkout fix_atgeometry command must be run from your at_geom_25Oct subdirectory. And if the directory already existed before your git clone command, you have to call git pull to pull the last version

lfarv · 2023-10-25T15:23:15Z

@simoneliuzzo: Python test

I reproduced exactly your commands, and I get no problem (see record below). On which platform did you run?
If it still fails, you can try as you did for Matlab: clone, change branch, pull and build.

bijou:devtest $ python3.9 -m venv atgeo
bijou:devtest $ source atgeo/bin/activate
(atgeo) bijou:devtest $ pip install --upgrade pip
Requirement already satisfied: pip in ./atgeo/lib/python3.9/site-packages (23.2.1)
Collecting pip
  Obtaining dependency information for pip from https://files.pythonhosted.org/packages/47/6a/453160888fab7c6a432a6e25f8afe6256d0d9f2cbd25971021da6491d899/pip-23.3.1-py3-none-any.whl.metadata
  Using cached pip-23.3.1-py3-none-any.whl.metadata (3.5 kB)
Using cached pip-23.3.1-py3-none-any.whl (2.1 MB)
Installing collected packages: pip
  Attempting uninstall: pip
    Found existing installation: pip 23.2.1
    Uninstalling pip-23.2.1:
      Successfully uninstalled pip-23.2.1
Successfully installed pip-23.3.1
(atgeo) bijou:devtest $ pip install git+https://github.com/atcollab/at.git@fix_atgeometry
Collecting git+https://github.com/atcollab/at.git@fix_atgeometry
  Cloning https://github.com/atcollab/at.git (to revision fix_atgeometry) to /private/var/folders/kc/97y2h1wj6sl5gp6cbqzgmyfc000105/T/pip-req-build-374vtx89
  Running command git clone --filter=blob:none --quiet https://github.com/atcollab/at.git /private/var/folders/kc/97y2h1wj6sl5gp6cbqzgmyfc000105/T/pip-req-build-374vtx89
  Running command git checkout -b fix_atgeometry --track origin/fix_atgeometry
  Switched to a new branch 'fix_atgeometry'
  branch 'fix_atgeometry' set up to track 'origin/fix_atgeometry'.
  Resolved https://github.com/atcollab/at.git to commit d843aa3b7a2d5646c83d5da4afded65b9aef28a3
  Installing build dependencies ... done
  Getting requirements to build wheel ... done
  Preparing metadata (pyproject.toml) ... done
Collecting numpy>=1.16.6 (from accelerator-toolbox==0.4.1.dev32+gd843aa3)
  Using cached numpy-1.26.1-cp39-cp39-macosx_10_9_x86_64.whl.metadata (61 kB)
Collecting scipy>=1.4.0 (from accelerator-toolbox==0.4.1.dev32+gd843aa3)
  Using cached scipy-1.11.3-cp39-cp39-macosx_10_9_x86_64.whl.metadata (60 kB)
Using cached numpy-1.26.1-cp39-cp39-macosx_10_9_x86_64.whl (20.6 MB)
Using cached scipy-1.11.3-cp39-cp39-macosx_10_9_x86_64.whl (37.3 MB)
Building wheels for collected packages: accelerator-toolbox
  Building wheel for accelerator-toolbox (pyproject.toml) ... done
  Created wheel for accelerator-toolbox: filename=accelerator_toolbox-0.4.1.dev32+gd843aa3-cp39-cp39-macosx_13_0_x86_64.whl size=1853834 sha256=ed4973c9c7434aca79975290785d0bd865d882277e30c46cb1215eb27e1a1e24
  Stored in directory: /private/var/folders/kc/97y2h1wj6sl5gp6cbqzgmyfc000105/T/pip-ephem-wheel-cache-i1epn219/wheels/a5/25/4b/a191292f1c614d4ab0955010de94974e161ba85ac34c6d43a8
Successfully built accelerator-toolbox
Installing collected packages: numpy, scipy, accelerator-toolbox
Successfully installed accelerator-toolbox-0.4.1.dev32+gd843aa3 numpy-1.26.1 scipy-1.11.3
(atgeo) bijou:devtest $ python
Python 3.9.18 (main, Aug 24 2023, 21:32:48) 
[Clang 14.0.3 (clang-1403.0.22.14.1)] on darwin
Type "help", "copyright", "credits" or "license" for more information.
>>> import at
>>> 
(atgeo) bijou:devtest $

lfarv · 2023-10-25T15:31:19Z

@simoneliuzzo : I may have an explanation for you python problem: in the command

/usr/local/bin/python3.9 /Users/liuzzo/PycharmProjects/testatgeometry/testgeom.py

you are not calling the python from your virtual environment, but the system one. AT is not installed there (or badly, it seems: File "/Users/liuzzo/beamdyn/at/pyat/at/__init__.py", line 3 what is this?).

While your virtual environment is activated, just call

python /Users/liuzzo/PycharmProjects/testatgeometry/testgeom.py

simoneliuzzo · 2023-10-25T16:29:21Z

Dear @lfarv ,

python solved. I forgot to set the interpreter.

now it is ok.

simoneliuzzo · 2023-10-25T16:30:59Z

contrary to atgetgeometry, python get_geometry, does not accept refpts as input. it is usefull, sometimes, to get the coordinates of a single point or a few, along the lattice.

simoneliuzzo · 2023-10-25T16:38:27Z

I did a gitshow in my at_geom_25Oct (named so, since I cloned it fresh just for this test.

! git show

commit d843aa3 (HEAD -> fix_atgeometry, origin/fix_atgeometry)�[m

Author: Laurent Farvacque laurent.farvacque@esrf.fr�[m

Date: Tue Oct 24 17:44:21 2023 +0200

Renamed _EPSIL to _GEOMETRY_EPSIL

I clean my path before doing the test. And add only the at in at_geom_25Oct

I check edit atgeometry and it seems the one of this branch.

simoneliuzzo · 2023-10-25T16:39:37Z

The lattice I use is may be an old one.

/machfs/liuzzo/EBS/beamdyn/matlab/optics/sr/S28F_restart/betamodel.mat

lfarv · 2023-10-25T18:45:24Z

@simoneliuzzo: even old lattices should work! I corrected the θ threshold detecting a full ring.
It now works with you lattice.
Thnaks for your tests, back to you.

simoneliuzzo · 2023-10-25T20:09:28Z

Thanks @lfarv ,

also matlab is ok for centered option.

lfarv · 2023-10-26T09:55:56Z

contrary to atgetgeometry, python get_geometry, does not accept refpts as input. it is usefull, sometimes, to get the coordinates of a single point or a few, along the lattice.

That's done. @simoneliuzzo, could you approve again? Thanks

simoneliuzzo · 2023-10-26T10:12:33Z

Checked. ok.
Blue is centered all elements
orange not centered at quads
green not centered at BPMs

lfarv added Matlab For Matlab/Octave AT code bug fix labels Oct 5, 2023

lfarv requested a review from simoneliuzzo October 5, 2023 09:42

lfarv added 2 commits October 5, 2023 11:45

correct the 'centered' flag

c941f13

Update the version number

f457fb7

lfarv force-pushed the fix_atgeometry branch from ddcbd8a to f457fb7 Compare October 5, 2023 09:46

lcarver reviewed Oct 9, 2023

View reviewed changes

lfarv added 2 commits October 23, 2023 16:46

Get center for full rings, half rings and single cells

524140a

Small fixes in type hints and docstrings

23af266

lfarv added the Python For python AT code label Oct 23, 2023

Renamed _EPSIL to _GEOMETRY_EPSIL

d843aa3

lfarv requested a review from swhite2401 October 24, 2023 15:47

corrected epsilon in atgeometry

1206d21

simoneliuzzo approved these changes Oct 25, 2023

View reviewed changes

added a refpts keyword in get_geometry

ef922a6

simoneliuzzo approved these changes Oct 26, 2023

View reviewed changes

lfarv merged commit 52298a9 into master Oct 26, 2023
31 checks passed

lfarv deleted the fix_atgeometry branch October 26, 2023 11:18

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Bug fix in atgeometry #667

Bug fix in atgeometry #667

lfarv commented Oct 5, 2023

simoneliuzzo commented Oct 9, 2023

lcarver Oct 9, 2023

simoneliuzzo Oct 9, 2023

simoneliuzzo commented Oct 10, 2023

lfarv commented Oct 10, 2023 •

edited

simoneliuzzo commented Oct 10, 2023

lfarv commented Oct 10, 2023

swhite2401 commented Oct 11, 2023

carmignani commented Oct 11, 2023

lfarv commented Oct 11, 2023

lfarv commented Oct 23, 2023

swhite2401 commented Oct 24, 2023

lfarv commented Oct 24, 2023

simoneliuzzo commented Oct 25, 2023

simoneliuzzo commented Oct 25, 2023

lfarv commented Oct 25, 2023

lfarv commented Oct 25, 2023

lfarv commented Oct 25, 2023 •

edited

simoneliuzzo commented Oct 25, 2023

simoneliuzzo commented Oct 25, 2023

simoneliuzzo commented Oct 25, 2023

simoneliuzzo commented Oct 25, 2023

lfarv commented Oct 25, 2023

simoneliuzzo commented Oct 25, 2023

lfarv commented Oct 26, 2023

simoneliuzzo commented Oct 26, 2023

Bug fix in atgeometry #667

Bug fix in atgeometry #667

Conversation

lfarv commented Oct 5, 2023

simoneliuzzo commented Oct 9, 2023

lcarver Oct 9, 2023

Choose a reason for hiding this comment

simoneliuzzo Oct 9, 2023

Choose a reason for hiding this comment

simoneliuzzo commented Oct 10, 2023

lfarv commented Oct 10, 2023 • edited

simoneliuzzo commented Oct 10, 2023

lfarv commented Oct 10, 2023

swhite2401 commented Oct 11, 2023

carmignani commented Oct 11, 2023

lfarv commented Oct 11, 2023

lfarv commented Oct 23, 2023

swhite2401 commented Oct 24, 2023

lfarv commented Oct 24, 2023

simoneliuzzo commented Oct 25, 2023

simoneliuzzo commented Oct 25, 2023

lfarv commented Oct 25, 2023

lfarv commented Oct 25, 2023

lfarv commented Oct 25, 2023 • edited

simoneliuzzo commented Oct 25, 2023

simoneliuzzo commented Oct 25, 2023

simoneliuzzo commented Oct 25, 2023

simoneliuzzo commented Oct 25, 2023

lfarv commented Oct 25, 2023

simoneliuzzo commented Oct 25, 2023

lfarv commented Oct 26, 2023

simoneliuzzo commented Oct 26, 2023

lfarv commented Oct 10, 2023 •

edited

lfarv commented Oct 25, 2023 •

edited