Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Adjusted approach to phase retrieval #36

Open
jacione opened this issue Oct 24, 2024 · 8 comments
Open

Adjusted approach to phase retrieval #36

jacione opened this issue Oct 24, 2024 · 8 comments

Comments

@jacione
Copy link
Contributor

jacione commented Oct 24, 2024

When the various algorithms are defined in op_flow.py, the propagations are considered as separate from the algorithms themselves:

algs = {'ER': ('to_reciprocal_space', 'modulus', 'to_direct_space', 'er'),
'HIO': ('to_reciprocal_space', 'modulus', 'to_direct_space', 'hio'),
'ERpc': ('to_reciprocal_space', 'pc_modulus', 'to_direct_space', 'er'),
'HIOpc': ('to_reciprocal_space', 'pc_modulus', 'to_direct_space', 'hio'),
}

I wonder if it might be easier to implement other algorithms if that were not the case. I'm specifically thinking of Table 1 in this paper, in which the entire update function is defined using projection operators.

@jacione jacione closed this as not planned Won't fix, can't repro, duplicate, stale Oct 24, 2024
@jacione jacione reopened this Oct 24, 2024
@brussel13
Copy link
Contributor

The comment above the alg definition says "four steps". Could an algorithm be more the 4 steps?

@jacione
Copy link
Contributor Author

jacione commented Oct 24, 2024

The steps are basically just "project, constrain, project, constrain", which I think is pretty foundational to any kind of alternating projection algorithm. However, some algorithms have nuances that allow communication between the "states" at various points in that cycle, like how HIO combines the results of the object before and after the modulus constraints.

What I'm wondering is if it might be better to think of the Fourier constraint as a single operation consisting of a forward propagation, modulus replacement, and back propagation. Doing so would allow us to define any update function as a linear combination of three constraint operators (identity, modulus, support)

@bfrosik
Copy link
Contributor

bfrosik commented Oct 24, 2024

Refer to (https://cohere-dev.readthedocs.io/en/latest/for_developers.html#adding-new-algorithm)
The documentation explains step by step how to add a new algorithm.
We can add the mentioned algorithms (I think we had SF at some point) any time.
Let's get input from Ross/Wonsuk about it.

The four steps are the often referred to in papers regarding BCDI, as fft, modulus projector, pfft, and ER/HIO (or other).

@brussel13
Copy link
Contributor

brussel13 commented Oct 24, 2024 via email

@bfrosik
Copy link
Contributor

bfrosik commented Oct 25, 2024

It doesn't have to be four steps. Any number of functions that are associated with the key (for example ER) will run.

@jacione
Copy link
Contributor Author

jacione commented Oct 25, 2024

After giving the matter some more thought, I think I may have misunderstood my own concern. I would reformulate it in this way:

I wonder if putting some separation between the constraints and the Rec/CoupledRec class might make it easier to implement new algorithms.

The main thing that I'm looking at is the difference map algorithm, which involves applying the RS modulus constraint twice, one of which is on an array that is not self.ds_image. Certainly that's not impossible to implement with some extra functions, but I'm just wondering if there might be a more elegant way. As a rough idea maybe something like this?

def to_reciprocal_space(self, arr=None):
    if arr is None:
        arr = self.ds_image
    self.rs_amplitudes = devlib.ifft(arr)

I don't feel super strongly about this, but it might be worth thinking about.

@brussel13
Copy link
Contributor

brussel13 commented Oct 25, 2024 via email

@jacione
Copy link
Contributor Author

jacione commented Oct 29, 2024

This makes sense. I was playing with it last week and came up with the following constraints for each of the algorithms listed in that Marchesini paper, but when I actually put them into cohere it didn't seem to work. They were just applying the initial support constraint and nothing else (not even modulus), so I don't know if there's just something weird going on in op_flow. I'm currently troubleshooting other things, but these might be a good starting point if you do end up working on it before I get back to it.

# Error reduction (already implemented)
def er(self):
    self.ds_image = self.ds_image_proj * self.support

# Solvent flipping (currently implemented as "new_alg")
def sf(self):
    self.ds_image = 2.0 * (self.ds_image_proj * self.support) - self.ds_image_proj

# Hybrid input-output (already implemented)
def hio(self):
    combined_image = self.ds_image - self.ds_image_proj * self.params['hio_beta']
    self.ds_image = devlib.where((self.support > 0), self.ds_image_proj, combined_image)

# Difference map
def dm(self):
    gs = -1/self.params['hio_beta']
    gm = 1/self.params['hio_beta']
    part_s = self.support * ((1+gs)*self.ds_image_proj - gs*self.ds_image)
    part_m_pre = (1+gm)*self.ds_image*self.support - gm*self.ds_image
    part_m = devlib.ifft(self.iter_data * part_m_pre / devlib.absolute(part_m_pre))
    self.ds_image = self.ds_image + self.params['hio_beta'] * (part_s - part_m)

# Averaged successive reflections
def asr(self):
    self.ds_image = self.ds_image + (2*self.ds_image_proj - self.ds_image)*self.support - self.ds_image_proj

# Hybrid projection reflection
def hpr(self):
    self.ds_image = 0.5*(self.support*((1+self.params['hio_beta'])*self.ds_image_proj - self.ds_image) + self.ds_image + (1-self.params['hio_beta'])*self.ds_image_proj

# Relaxed averaged alternating reflectors
def raar(self):
    self.ds_image = self.ds_image_proj + self.params['hio_beta']* ((2*self.ds_image_proj - self.ds_image)*self.support - 2*self.ds_image_proj + self.ds_image)

Obviously, it won't be the hio_beta parameter for all of these, but many of them have their own weighting parameter.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

3 participants