Adjusted approach to phase retrieval

[jacione]

When the various algorithms are defined in op_flow.py, the propagations are considered as separate from the algorithms themselves: https://github.com/AdvancedPhotonSource/cohere/blob/7367cd151004253b398b464d349171cff1e43948/cohere_core/controller/op_flow.py#L7-L11

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.

[brussel13]

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

[jacione]

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]

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]

I was asking if it can be more or less than four steps. separating them is important. one might not use the fft to do the propagation but the modulus could still be the same.

---

From: Barbara Frosik @.> Sent: Thursday, October 24, 2024 11:43 PM To: AdvancedPhotonSource/cohere @.> Cc: Harder, Ross J. @.>; Comment @.> Subject: Re: [AdvancedPhotonSource/cohere] Adjusted approach to phase retrieval (Issue #36)

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.  ZjQcmQRYFpfptBannerStart This Message Is From an External Sender This message came from outside your organization.

ZjQcmQRYFpfptBannerEnd

Refer to (https://cohere-dev.readthedocs.io/en/latest/for_developers.html#adding-new-algorithmhttps://urldefense.us/v3/__https://cohere-dev.readthedocs.io/en/latest/for_developers.html*adding-new-algorithm__;Iw!!G_uCfscf7eWS!fcyJZKGeNQ5abjmezK9RVZ3M2J71ZYaGsZtyAyx7nej8gLD8kg0b8FH5qwEVFyu_vPjER9NzcclzGkqLZm77xHNCcw$) 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).

— Reply to this email directly, view it on GitHubhttps://urldefense.us/v3/__https://github.com/AdvancedPhotonSource/cohere/issues/36*issuecomment-2436387978__;Iw!!G_uCfscf7eWS!fcyJZKGeNQ5abjmezK9RVZ3M2J71ZYaGsZtyAyx7nej8gLD8kg0b8FH5qwEVFyu_vPjER9NzcclzGkqLZm7QtAPgsA$, or unsubscribehttps://urldefense.us/v3/__https://github.com/notifications/unsubscribe-auth/ABZ6SOBKVQQONQISBALRPD3Z5FSYVAVCNFSM6AAAAABQR5P7CSVHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMZDIMZWGM4DOOJXHA__;!!G_uCfscf7eWS!fcyJZKGeNQ5abjmezK9RVZ3M2J71ZYaGsZtyAyx7nej8gLD8kg0b8FH5qwEVFyu_vPjER9NzcclzGkqLZm5KlqAd4g$. You are receiving this because you commented.Message ID: @.***>

[bfrosik]

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

[jacione]

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]

Nick,

We had some discussion about this. At least for now, we decided that we wanted to keep the algorithm operators argument free. It makes the algorithm definition in op_flow clear. If there were arguments possible then we would need some way of including those in the algorithm definition and it wasn't completely obvious how to do that. We can start working on an DM implementation. Just need another "to_reciprocal_space" operator in the class to propagate the other direct space image used in the algorithm.

---

From: Nick Porter @.> Sent: Friday, October 25, 2024 7:46 PM To: AdvancedPhotonSource/cohere @.> Cc: Harder, Ross J. @.>; Comment @.> Subject: Re: [AdvancedPhotonSource/cohere] Adjusted approach to phase retrieval (Issue #36)

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 ZjQcmQRYFpfptBannerStart This Message Is From an External Sender This message came from outside your organization.

ZjQcmQRYFpfptBannerEnd

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.

— Reply to this email directly, view it on GitHubhttps://urldefense.us/v3/__https://github.com/AdvancedPhotonSource/cohere/issues/36*issuecomment-2438435275__;Iw!!G_uCfscf7eWS!ZB5gk70IOyFldPhTYicE6YhIYGfc77WQ81sSHYnj-T6DaqKQWQ0Jf7v-poxiOVL9oH5qX0KuJ3if5jIZjllFOeagQQ$, or unsubscribehttps://urldefense.us/v3/__https://github.com/notifications/unsubscribe-auth/ABZ6SOFRE37DUR444JTHKJLZ5J7VVAVCNFSM6AAAAABQR5P7CSVHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMZDIMZYGQZTKMRXGU__;!!G_uCfscf7eWS!ZB5gk70IOyFldPhTYicE6YhIYGfc77WQ81sSHYnj-T6DaqKQWQ0Jf7v-poxiOVL9oH5qX0KuJ3if5jIZjlmYZK_MAw$. You are receiving this because you commented.Message ID: @.***>

[jacione]

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.