atcollab / at

Accelerator Toolbox
Apache License 2.0
48 stars 31 forks source link

add EApertures and RApertures #635

Closed oscarxblanco closed 11 months ago

oscarxblanco commented 11 months ago

Dear all,

this Pull Request solves an issue on the aperture name case for the first element in a ring when saving a .mat file

For example, test.mat was saved and read back; the title case sets Eapertures to a lower case 'a' .

>>> r1 = at.load_mat('test.mat',keep_all=True)
>>> r1[0]
Marker('ringtest', Beam_Current=0.0, Eapertures=array([0.0075, 0.005 ]), ...

In this pull request the title case dictionary is modified. Also, EApertures and RApertures are included in _param_to_lattice.

oscarxblanco commented 11 months ago

Dear @swhite2401 , is this modification ok ? who should review it ?

swhite2401 commented 11 months ago

Hello @oscarxblanco, I could not reproduce the problem with the following code:

import at
ring=at.Lattice([at.Marker('test', Beam_Current=0.0, EApertures=[0.0075, 0.005 ])], energy=6e9, periodicity=1)
at.save_lattice(ring,'test.mat')
ring2 = at.load_mat('test.mat')
print(ring2[0])

Could send a full example?

oscarxblanco commented 11 months ago

@swhite2401 It is an issue that appears when the flag keep_all=True is used.

There are several scenarios, the first is when the lattice is created from a python script, then written, and then read into python. A second scenario is when the lattice is already in .mat format with aperture fields on the first element. A third scenario appears when the lattice is loaded in matlab and then saved into another .mat file.

I'm preparing the test for these scenarios and any other if I think I missed any. In the meanwhile here is the output I get if I add the flag keep_all=True to your code :

import at
ring=at.Lattice([at.Marker('test',
                            Beam_Current=0.0,
                            EApertures=[0.0075, 0.005 ])],
                 energy=6e9,
                 periodicity=1
                 )
print(ring[0])
at.save_lattice(ring,'test.mat')
ring2 = at.load_mat('test.mat',keep_all=True)
print(ring2[0])

OUTPUT:

Marker:
        FamName : test
        Length : 0.0
        PassMethod : IdentityPass
        Beam_Current : 0.0
        EApertures : [0.0075 0.005 ]
Marker:
        FamName :
        Length : 0.0
        PassMethod : IdentityPass
        Particle : relativistic
        Beam_Current : 0.0
        Nbunch : 1.0
        Periodicity : 1
swhite2401 commented 11 months ago

Isn't this the normal behavior -> RingParam is kept in a marker and EAperture has the correct case right? What is strange is that the RingParam Marker has no FamName... I think this can be a problem

oscarxblanco commented 11 months ago

@swhite2401 In the second print EApertures does not longer appear.

And yes, FamName should not disappear either.

lfarv commented 11 months ago

@oscarxblanco : you are not correcting the real problem. The problem is the handling RingParam elements, not a problem of apertures.

I'll look at that.

lfarv commented 11 months ago

What is strange is that the RingParam Marker has no FamName... I think this can be a problem

That's normal: the lattice has not been given any name, so it's an empty string

lfarv commented 11 months ago

@oscarxblanco and @swhite2401:

For me, everything looks normal:

ring=at.Lattice([at.Marker('test',
                            Beam_Current=0.0,
                            EApertures=[0.0075, 0.005 ])],
                 energy=6e9,
                 periodicity=1
                 )
at.save_lattice(ring,'test.mat')
ring2 = at.load_mat('test.mat',keep_all=True)
print(ring2[0])
print(ring2[1])

OUTPUT:

Marker:
    FamName : 
    Length : 0.0
    PassMethod : IdentityPass
    Particle : relativistic
    Beam_Current : 0.0
    Nbunch : 1.0
    Periodicity : 1
Marker:
    FamName : test
    Length : 0.0
    PassMethod : IdentityPass
    Beam_Current : 0.0
    EApertures : [0.0075 0.005 ]

The first element is the conversion of the RingParam pseudo-element into a Marker containing the lattice parameters, the 2nd element is the correct marker with the right EAperture field. And the extra Beam_Current attribute you assigned it.

Where do you see a problem?

I recall: By default (keep_all=False), PyAT takes lattice attributes from the Matlab RingParam elements and discards them since they are useless. With keep_all=True, PyAT keeps the RingParam elements as Markers. This is useful to keep the number of elements identical when comparing Matlab and python results.

oscarxblanco commented 11 months ago

Dear @lfarv , I partially disagree with your point of view of the keep_all=True usage. Many users still prefer matlab and some are trying out pyat, therefore, keep_all=True is not only limited to comparison of results or cross check, it is the way to give them the chance to work in one or both environments and at the same time keeping identical indexing. Having the same indexing and lattice length allows to reduce the chances of errors when sharing results.

The fact that RingParam is treated differently in matlab AT and pyat produces differences in .mat files that are a common way to communicate both systems. One could then find oneself in the situation where identical lattices saved in different environment will differ and therefore it requires that one keeps track of the environment used to produced the results or to write additional custom code so that the differences in the first element are corrected. For example in this particular case, if I read your 'test.mat' file in my matlab AT session I will have to make sure to create again the EApertures you already wrote in your python script for your simulations, and correct any other modification when you saved in pyat and I load in matlab AT.

## python
...
ring=at.Lattice([at.Marker('test',
                            Beam_Current=0.0,
                            EApertures=[0.0075, 0.005 ])], # your EApertures in pyat !
                 energy=6e9,
                 periodicity=1
                 )
...
>> %% matlab
>> load('test.mat')
>> RING
RING =
  2×1 cell array
    {1×1 struct}
    {1×1 struct}
>> RING{1} % my first element does not longer have EApertures :\
ans = 
  struct with fields:
         FamName: ''
          Length: 0
      PassMethod: 'IdentityPass'
        Particle: [1×1 struct]
     Periodicity: 1
    Beam_Current: 0
          Nbunch: 1
          Energy: 6.0000e+09
           Class: 'RingParam'
>> RING{2}
ans = 
  struct with fields:
         FamName: 'test'
          Length: 0
      PassMethod: 'IdentityPass'
    Beam_Current: 0
      EApertures: [0.0075 0.0050]
           Class: 'Marker'

However, I understand now that the prefered way to work in pyat is to remove RingParam from the sequence, and therefore the first element I show won't exist. Matlab users should never add anything to RingParam that is not supported by pyat (in this case EApertures), and about the indexing, we would always have one less element in pyat wrt matlab AT.

Developing a use example a little bit further, when saving results of let's say momentum acceptance where I need to keep a list of indexes and two arrays (positive and negative local momentum acceptance) with the same length of the indexes, I would need also to know if the saved indexing was produced in matlab AT or pyat, so that I substract 1 from matlab to python, and another one from the default keep_all false. In that way anyone using the touschek_lifetime routine on matlab or python would have the same result.

Thank you for your input @lfarv and @swhite2401 , if there is no more comments I will close this pull request without merging by the end of the week.

swhite2401 commented 11 months ago

I would like to insist on the fact that having the RingParam element with empty FamName is quite misleading, I think it would help a lot comprehension if it was named RingParam for lattice without name instead. In fact why not always call it RingParam, is there any reason for that?

lnadolski commented 11 months ago

I would like to insist on the fact that having the RingParam element with empty FamName is quite misleading, I think it would help a lot comprehension if it was named RingParam for lattice without name instead. In fact why not always call it RingParam, is there any reason for that?

I fully agree with Simon. An empty name is not a good practice and can also be a source of failure for some functions.

lfarv commented 11 months ago

@swhite2401 and @lnadolski: Sorry, but the FamName of the RingParam element is by convention used to store the lattice name. This name is used in plot titles, and possibly other places. In PyAt, it is set as the Lattice name attribute. In both Matlab and python, the default is an empty string. If you don't like empty strings, I'm ready to set anything else for the default lattice name (to be set in the python Lattice constructor and in the Matlab RingParam constructor). Any proposals? (I think a plot title like "ringParam" is not a good idea!)

swhite2401 commented 11 months ago

Can't you simply add and attribute LatticeName to store the lattice name and always set FamName=RingParam?

lfarv commented 11 months ago

Can't you simply add and attribute LatticeName to store the lattice name and always set FamName=RingParam?

It could have been done that way… But:

But again, if you think of a better default name for lattices, it's easy to change, without any compatibility problem: lattices created and stored with an empty name will keep it.

lfarv commented 11 months ago

For @oscarxblanco, some explanation on the RingParam element: In Matlab, there is no lattice object, instead, all the attributes of a lattice are stored as attributes of a RingParam pseudo-element. Practically, all the RingParam attributes are mapped in python to Lattice attributes. On the other side, when creating a .mat file, all the python Lattice attributes are set as attributes of a RingParam element prepended to the lattice elements. There are a few predefined attributes (Energy, Particle…), but the mapping is also done for any extra user-defined attribute. So:

Your example shows exactly as it works: you lattice has a single element: a Marker with an EAperture and a Beam_current attribute (???). and you find it exactly identical (same name, same attributes) in Matlab. It's in the 2nd position in Matlab because of the presence of the RingParam added by python in the .mat file and necessary in Matlab, but it's strictly identical. And it will be also be in 2nd position when read again in python with keep_all=True.

swhite2401 commented 11 months ago

@lfarv this PR was opened because of a problem and the minimal solution you propose is not satisfactory in my opinion: the source of the confusion is still present, undocumented and we will face similar questions in the future.

I take note that you do not want to change your design even though this was introduced recently and I strongly doubt anyone is making use of the lattice name feature: AT users generally set the .mat file name to be the lattice name.

I can propose several actions to sort this out:

  1. Document this feature: this is absolutely necessary, and not just in the lattice_object docstring. A full manual on the lattice object is really needed for everyone to follow what you have implemented without having to spend hours looking at the source
  2. Set the FamName of the ringparam element to the file name without extension by default
  3. Introduce a RingParam element in pyAT and set the Class attribute to RingParam, such that this element is clearly identified as what it is
lfarv commented 11 months ago

@oscarxblanco : Sorry, I did not say I don't want to change anything, but just said it worked as foreseen.

this PR was opened because of a problem

I still do not see what is this problem ! In all examples you showed, you create a lattice with a single element, and then whatever you do, saving a file, reading it in Matlab or python, you single element is perfectly restored, strictly identical. Only in Matlab, its prepended with a RingParam element which is automatically generated because Matlab needs it. You must not care about its contents, it's private Matlab stuff, and this element is fully transparent in tracking.

To better understand, could you say exactly:

  1. what looks wrong to you,
  2. what you would expect

this was introduced recently and I strongly doubt anyone is making use of the lattice name feature

I looked in history, and found:

So I can assure that there is nothing new there ! The name of a lattice is different from a file name: it can be used to identify specific lattice features, versioning, or else. Whether you use it or not, it's up to you, it's completely optional.

In the mean time, I'll look at you suggestions.

Ideally, this stream should better be in 'issues' or 'discussions', but I do not know if it can be moved…

swhite2401 commented 11 months ago

Ok in fact it is not new I was wrong, I got sidetracked by the fact that it suddenly became mandatory for other functions to work. Whether people use it or not is another problem, but honestly I don't think they are and changing it will have little or no consequences.

The algorithm is perfectly fine (at least for me) so in this sense there is not problem. My concern is that this ringparam element is not documented (or is it??) or self descriptive and users are getting confused by it which I can perfectly understand! So this is why I proposed these changes that will not affect any of its functionalities but should at least provide all the information to identify this element and understands what it is used for.

I agree this should be moved to issues....but I have no idea on how to do this....

lfarv commented 11 months ago

I got sidetracked by the fact that it suddenly became mandatory

Again, I'm lost: nothing became mandatory: a RingParam element is not. It is used to improve performance by caching some values (energy for instance) which need a full scan of the lattice to be retrieved, or to store some parameters (name, particle) which all have default values in case it's missing. Missing a RingParam just lowers slightly the performance and prevents using recent or future features (non-relativistic tracking).

And nothing changed (and not suddenly) since long: the RingParam element is stored by python in .mat files since the very beginning (~2018). In Matlab, the last modification involving RingParam is the addition of atSetRingProperties/atGetRingProperties in november 2021.

What is the sudden change you are referring to?

oscarxblanco commented 11 months ago

Dear @lfarv , the explanation you gave is clear. RingParam is not an element.

Dear all, by tomorrow I will close this Pull Request without merging; would you need or prefer to move this to a discussion ? If so, how should I do it ? I won't like to cut out the discussion if you think it is important.

Best regards, o

swhite2401 commented 11 months ago

Well you told us that FamName of the RingParam element could not be anything (not RingParam as suggested) because it was used for plotting.... isn't this new? It seems we are just getting really confused on very simple issues...

I would like to suggest that we focus practical issues: how do you we make it clearer to users what this RingParam element is? what it is used for? and how do we clearly identify it as an non standard element?

I proposed something 3 mitigation. Are these relevant and can we implement them? Do you have other proposals? I have no strong opinion I am just trying to help here I you think it is perfect like this and nothing should be done that is fine by me but I insist that documentation on this point is required.