@cgohlke I started this PR so we can discuss the components module. I have been discussing with Leonel about what method is the best for calculating the fraction of two components and we arrived to the conclusion that both the more "mathematical" solution (solving system of equations) and the more "graphical" approach (moving cursor and creating histogram) should both be available. I started implementing the mathematical approach until the cursors module is ready.

Regarding this approach, I tried both with solving linear equations and the alternative by projecting points to a line between components and then calculating euclidean distance, which is inversely proportional to the fraction of the component. I tested and this last method was fastest than solving a system of equations, so I went for that approach, at least for now so we can start testing this implementations. The approach for solving equations is also necessary, particularly for calculating fractions of three and four components, so it will also hav to be implemented, but maybe if it is faster we can keep this implementation for two components and use the other for more components.

Let me know what do you think and we can continue this discussion since I have other issues I would like to ask you related to this PR.

This looks good. I left some comments.

"graphical" approach (moving cursor and creating histogram)

Can you explain this.

this last method was fastest than solving a system of equations

If you have a large number of phasor coordinates and are OK with a limited precision, you could calculate a 2D lookup table of fractions with reduced resolution, and then index that table with the digitized phasor coordinates. The lookup table can be calculated using your two_fractions_from_phasor function...

@cgohlke thanks for all the comments, they helped a lot! I think I addressed all of them, and some of them I replied in the same thread as to make it easier to follow.

"graphical" approach (moving cursor and creating histogram)

This is a method that Leonel explained to me, that he believes is what SimFCS does. This method is also described in this paper:
DOI: 10.1021/acs.jpca.9b07880

Briefly, the idea is that by moving the cursor along a line in discrete steps and then creating a histogram along the line in each step helps to take into account the dispersion in all directions, instead of the projection to the line that only considers one direction. Maybe we can discuss this next Monday with Leonel.

I also like the idea of a 2D lookup table. Maybe it's another thing to discuss on monday if it makes sense to compromise a little of resolution for speed, at least to have that option I think would be great.

If I may suggest some minor cleanup of the tutorial, mainly separating variable definitions from calculations and plotting by empty lines. Also, how about combining the two subplots since the histograms do not overlap much?

diff --git a/tutorials/phasorpy_components.py b/tutorials/phasorpy_components.py
index f6154dc..e2070d9 100644
--- a/tutorials/phasorpy_components.py
+++ b/tutorials/phasorpy_components.py
@@ -27,12 +27,15 @@ from phasorpy.plot import PhasorPlot
 
 frequency = 80.0
 components_lifetimes = [8.0, 1.0]
+component_fractions = [0.25, 0.75]
+
 real, imag = phasor_from_lifetime(
-    frequency, components_lifetimes, [0.25, 0.75]
+    frequency, components_lifetimes, component_fractions
 )
 components_real, components_imag = phasor_from_lifetime(
     frequency, components_lifetimes
 )
+
 plot = PhasorPlot(frequency=frequency, title='Combination of two components')
 plot.plot(components_real, components_imag, fmt='o-')
 plot.plot(real, imag)
@@ -47,7 +50,8 @@ plot.show()
     fraction_of_first_component,
     fraction_of_second_component,
 ) = two_fractions_from_phasor(real, imag, components_real, components_imag)
-print(f'Fraction of first component: {fraction_of_first_component:.3f}')
+
+print(f'Fraction of first component:  {fraction_of_first_component:.3f}')
 print(f'Fraction of second component: {fraction_of_second_component:.3f}')
 
 # %%
@@ -60,9 +64,10 @@ print(f'Fraction of second component: {fraction_of_second_component:.3f}')
 real, imag = numpy.random.multivariate_normal(
     (0.6, 0.35), [[8e-3, 1e-3], [1e-3, 1e-3]], (100, 100)
 ).T
+
 plot = PhasorPlot(
     frequency=frequency,
-    title='Phasor with contibution of two known components',
+    title='Phasor with contribution of two known components',
 )
 plot.hist2d(real, imag, cmap='plasma')
 plot.plot(*phasor_from_lifetime(frequency, components_lifetimes), fmt='o-')
@@ -77,15 +82,26 @@ plot.show()
     fraction_from_first_component,
     fraction_from_second_component,
 ) = two_fractions_from_phasor(real, imag, components_real, components_imag)
-fig, ax = plt.subplots(2, 1, figsize=(8, 8))
-ax[0].hist(fraction_from_first_component.flatten(), range=(0, 1), bins=100)
-ax[0].set_title('Histogram of fractions of first component')
-ax[0].set_xlabel('Fraction of first component')
-ax[0].set_ylabel('Counts')
-ax[1].hist(fraction_from_second_component.flatten(), range=(0, 1), bins=100)
-ax[1].set_title('Histogram of fractions of second component')
-ax[1].set_xlabel('Fraction of second component')
-ax[1].set_ylabel('Counts')
+
+fig, ax = plt.subplots()
+ax.hist(
+    fraction_from_first_component.flatten(),
+    range=(0, 1),
+    bins=100,
+    alpha=0.75,
+    label='First',
+)
+ax.hist(
+    fraction_from_second_component.flatten(),
+    range=(0, 1),
+    bins=100,
+    alpha=0.75,
+    label='Second',
+)
+ax.set_title('Histograms of fractions of first and second component')
+ax.set_xlabel('Fraction')
+ax.set_ylabel('Counts')
+ax.legend()
 plt.tight_layout()
 plt.show()
 

Excelent! I like the modifications you did to the tutorial. I will update it with your suggestions.

Hi @cgohlke, I am having some trouble trying to make both the project_phasor_to_line and two_fractions_from_phasor to be able to handle inputs with several dimensions. I am not exactly sure if this will be necessary also. For example if the real and imag are of shape (512, 512, 4) because they have different harmonics or channels (4), then real_components and imag_components should be also of shape (2, 4) and the projection should be performed over the last dimension. I am trying but couldn't find a solution.

Do you think this is necessary? If so, do you have any suggestion of how to approach this with the current implementation?

I haven't looked in that detail yet. If you think supporting those cases will complicate the implementation or documentation too much, how about documenting the limitation in the docstring and adding a test that should ideally pass but currently fails (marked pytest.xfail).

Do you think this is necessary?

Let's leave it out of this PR. Most methods/algorithms in this module will not naturally support multiple channels and/or frequencies (no?). If we decide to re-implement this algorithm as a numpy ufunc in Cython, those cases will be supported automatically.

Totally agree, I will update documentation and tests and let you know when it is ready again for reviewing.

@cgohlke I think this is ready for merge, but please let me know if something else can be improved

For some reason I get the following doctest error:

        Examples
        --------
        >>> two_fractions_from_phasor(
        ...     [0.6, 0.5, 0.4], [0.4, 0.3, 0.2], [0.2, 0.9], [0.4, 0.3]
        ... ) # doctest: +NUMBER
        (array([0.44, 0.56, 0.68]), array([0.56, 0.44, 0.32]))

        """
        real_components = numpy.asarray(real_components)
        imag_components = numpy.asarray(imag_components)
        if real_components.shape != (2,):
>           raise ValueError(f'{real_components.shape=} != (2,)')
E           ValueError: real_components.shape=(2, 3) != (2,)

For some reason I get the following doctest error:

        Examples
        --------
        >>> two_fractions_from_phasor(
        ...     [0.6, 0.5, 0.4], [0.4, 0.3, 0.2], [0.2, 0.9], [0.4, 0.3]
        ... ) # doctest: +NUMBER
        (array([0.44, 0.56, 0.68]), array([0.56, 0.44, 0.32]))

        """
        real_components = numpy.asarray(real_components)
        imag_components = numpy.asarray(imag_components)
        if real_components.shape != (2,):
>           raise ValueError(f'{real_components.shape=} != (2,)')
E           ValueError: real_components.shape=(2, 3) != (2,)

That's strange, I can't seem to replicate that error. That comes from testing the example in the docstring? Could it be something related to the pytest.xfail test? I think is the only place I test with a real_components of shape (2,3)

Never mind. I misread the pytest output. It's the output of the known failure. All tests passed!