inconsistencies in algorithms modules

[pavlis]

I thought about posting this as an issues page, but I think it need a bit of online discussion as a preface to addressing the problem. The issue we have is that the standard mspass algorithms modules are not consistent in handling ensembles. Many of them do but many do not either. Worse, is that the documentation doesn't alway make it clear what data types an algorithm supports. The reason it is not a simple "fix this problem" is that some algorithms make sense only on certain kinds of data. e.g. rotation of a Seismogram object to R-T-Z from ENZ only makes sense really for Seismogram objects. (Obspy actually has a way to handle that be that is not the point.) Another example is that an ensemble stack algorithm is not at all the same thing as atomic sums done with reduce. Hence, a (hypothetical in this case since it doesn't currently exist in mspass) ensemble_stack function would only want to accept ensembles as inputs. Add to that the Gather version of ensembles and it gets even messier.

I think the things we need to do are:

1. Maybe we have one but I didn't immediately find it, but we need a decorator that could be used with atomic algorithms to simplify running them on ensembles. The only thing I see to do that is to adapt obspy algorithms where that kind of comes along for free. The issue is when the atomic algorithm is not obspy I've written a bunch of wrappers for individual algorithms and in the simplest cast the problem boils down to looping the algorithm over ensemble members. Done right may require a few details to be more robust.
2. We need to go through the docstrings of all the standard algorithms and make sure there is a :return: section in the docstring that clarifies what exactly the algorithm will return as well as what it accepts as inputs.
3. I'm planning a new user manual page with a table of algorithm function names, one line concept summary, and a link to the docstring. Seems that table would benefit from an addition column with a symbolic summary of inputs allowed and outputs emitted. By symbolic I mean the table contents would be strings like "T S T_e S_e" that would be translated to latex symbols. T=TimeSeries, S=Seismogram, T_e=TimeSeriesEnsemble, and S_e is Seismogram ensemble. Right now we don't have any direct support for the prototype Gather object but the table could define a standard symbol for that too

Please let me know what any of you think about this.

[wangyinz]

Yeah, I think they are all pretty good ideas. We should implement them. In item 1, I don't think we have to complicate it, at least for now. I think it is sufficient to just have a loop over members. The other more complex cases may not be that common and decorator will not be necessary in that case.

[pavlis]

Let's talk about this at our next meeting then. For item 1, I think you could write some generic wrapper that would have the processing function and the args to be passed to it. I don't exactly know how to do that. The extra fluff I would recommend is a sanity check on the type of arg0 - it needs to be a mspass ensemble or gather object.

[pavlis]

I ran across a subtle and obnoxious python anomaly on this issue that I think we have discussed previously. That is, there are two ways you can apply an atomic-level function to an ensemble e.

1. The fortran way with a subscript like operator:

for i in range(len(e.member)):
    e.member[i]=process(e.member[i])

2. The iterator way:

for d in e.member:
    d = process(d)

I don't know if it is a python anomaly or something about our pybind11 binding code, but they do not give the same answer. 1 works but 2 does nothing except consume cycles.

We should resolve this. The iterator way is natural to most python programmers and it produces mysterious results when it does nothing,.

[pavlis]

Working on this a bit more I remembered this is a weird python problem, not a pybind11 issue in our implementation. The iterator example I gave above illustrates the problem. If the loop is of that basic form all the loop does is continuously recycle the symbol d. I think one way to say it is that d has scope limited to the loop so python never associates it with the symbol e except to define the iterator sequence.

The really odd behavior comes when you use a loop like this one:

for d in e.member:
    d.rtoa()

Turns out that works because rota is a method of the object d points to so the action is pushed to the actual data object. Similar things work if you do constructs like:

for d in e.member:
    d['foo'] = 'bar'

In that case d actually worms its way into e to alter that data.

This is going to be very very confusing to almost anyone, Three actions:

1. We need to go through all our algorithm modules and make sure we don't have a hidden bug due to this oddity.
2. When this discussion turns into something concrete whomever does this needs to be very careful with this issue.
3. Since I'm in the middle of this documentation revision I need to find a place to document this subtle point and explain it.