Some phase tweaks... kind of work in progress.

Author: StephenThornquist

siffpy/core/siffreader.py

@@ -1946,6 +1946,7 @@ def register(
             kwargs.pop('nowarn')
         else:
             warnings.warn("\n\n \t Don't forget to fix the zplane alignment!!")
+
         if not self.opened:
             raise RuntimeError("No open .siff or .tiff")
 

siffpy/core/utils/circle_fcns.py

@@ -1,5 +1,6 @@
 # Functions for circularizing floats and ints
 from typing import Any, Tuple
+import warnings
 
 import numpy as np
 
@@ -180,6 +181,11 @@ def circ_corr(
     >>> 0.9131513426682732
     ```
     """
+    if (x.dtype == np.complex128) or (y.dtype == np.complex128):
+        warnings.warn(
+            "Passed complex-valued array to `circ_corr`. May be a mistake. \
+            Consider using `circ_corr_complex` instead."
+        )
     expd_1 = np.exp(1j*x)
     expd_2 = np.exp(1j*y)
 

siffpy/siffmath/flim/traces.py

@@ -193,7 +193,7 @@ def subtract_noise(
             presumed to be the same as the `FlimTrace`.
         """
 
-        # Subtracts a known amount of noise.
+        # Subtracts the noise if the value is already known.
         if n_photons is not None:
             if (
                 isinstance(n_photons, np.ndarray) 
@@ -229,6 +229,7 @@ def subtract_noise(
         # FLIMParams object.
         if self.FLIMParams is None:
             return
+        
         # Uses the value of the noise, assumes a uniform distribution,
         # and subtracts that proportion of photons from the intensity
         # and lifetime data.

siffpy/siffmath/phase/__init__.py

@@ -88,7 +88,7 @@ def phase_shift(x : np.ndarray, shift : PhaseTraceLike)->np.ndarray:
         idx = np.angle(np.exp(1j*shift[t]))*n_cols/(2*np.pi)
         whole = idx.astype(int)
         frac = idx - whole
-        shifted[:,t] = np.roll(np.abs(1-frac)*x[:,t], whole)
+        shifted[:,t] = np.roll((1-np.abs(frac))*x[:,t], whole)
         shifted[:,t] += np.roll(np.abs(frac)*x[:,t], int(whole+np.sign(frac)))
 
     return shifted

siffpy/siffmath/phase/correlate.py

@@ -0,0 +1,76 @@
+import numpy as np
+from typing import Iterable
+
+from siffpy.core.utils.circle_fcns import circ_corr_complex, running_circ_corr_complex
+from .traces import PhaseTrace
+
+def circular_correlate(x : PhaseTrace, y : PhaseTrace, method : str = 'Fisher', ) -> np.ndarray:
+    """
+    Circular correlation of two phase traces -- does NOT take error functions into
+    account. TODO: use the magnitude rather than just the angle for a more robust
+    estimate of error.
+
+    ## Arguments
+
+    x : PhaseTrace
+        The first phase trace to correlate.
+
+    y : PhaseTrace
+        The second phase trace to correlate.
+
+    method : str
+
+        The method to use for the correlation. Options are:
+        - 'Fisher'
+        - 'Jammalamadaka
+
+    ## Returns
+
+    np.ndarray
+        The circular correlation of the two phase traces. A single number.
+    """
+
+    return circ_corr_complex(
+        np.exp(1j*x.angle),
+        np.exp(1j*y.angle),
+        method = method
+    )
+
+def correlate(x : PhaseTrace, y : PhaseTrace, dts : Iterable, method : str = 'Fisher', mode : str = 'valid') -> np.ndarray:
+    """
+    Compute the circular cross-correlation of two phase traces, the correlation of `x` with a time-shifted
+    version of `y`. Does NOT take error functions or the magnitude of the phase into account.
+
+    ## Arguments
+
+    x : PhaseTrace
+        The first phase trace to correlate.
+
+    y : PhaseTrace
+        The second phase trace to correlate.
+
+    method : str
+
+        The method to use for the correlation. Options are:
+        - 'Fisher'
+        - 'Jammalamadaka
+
+    mode : str
+
+        The mode to use for the correlation. Options are:
+        - 'full'
+        - 'valid'
+        - 'same'
+
+    ## Returns
+
+    np.ndarray
+        The cross-circular-correlation of the two phase traces.
+    """
+    if not mode == 'valid':
+        raise NotImplementedError("Only 'valid' mode is currently supported")
+    
+    x_angle = np.exp(1j*x.angle)
+    y_angle = np.exp(1j*y.angle)
+
+    raise NotImplementedError("This function is not yet implemented")

siffpy/siffmath/phase/phase_estimates.py

@@ -5,7 +5,7 @@
 an argument vector_timeseries, which is a numpy array,
 and accept a keyword argument error_estimate, which is a boolean
 """
-from typing import Callable, Union, Optional, Tuple
+from typing import Callable, Union, Optional, Tuple, Any
 from enum import Enum
 
 import numpy as np
@@ -22,7 +22,8 @@
 )
 
 __all__ = [
-    'pva'
+    'pva',
+    'pva_flim',
 ]
 
 class PhaseErrorFunction(Enum):
@@ -101,6 +102,7 @@ def pva(
         time              : Optional[np.ndarray]              = None,
         error_function    : Optional[Union[Callable,str]]     = 'relative_magnitude',
         filter_fcn        : Optional[Union[Callable,str]]     = None,
+        angle_coords      : Optional[np.ndarray[Any, Any]]              = None,
         **kwargs
     ) -> PhaseTrace:
     """
@@ -136,6 +138,15 @@ def pva(
     filter_fcn : function
 
         A function to apply to the time series before computing the PVA.
+
+    angle_coords : np.ndarray
+
+        The angular coordinates of the first dimension of the vector timeseries,
+        i.e. the angle that each vector component corresponds to. If not provided,
+        will assume that the angle is linearly spaced between -pi and pi. Can
+        provide complex numbers on the unit circle or floats corresponding to an
+        angle in radians.
+    
     """
 
     if normalized:
@@ -144,8 +155,7 @@ def pva(
         max_val = sorted_vals[:,int(sorted_vals.shape[-1]*(1.0-1.0/20))]
         vector_timeseries = ((vector_timeseries.T - min_val)/(max_val - min_val)).T
 
-    angle_coords = np.exp(np.linspace(np.pi, -np.pi, vector_timeseries.shape[0])*1j) # it goes clockwise.
-    if isinstance(vector_timeseries, FluorescenceTrace):
+    if isinstance(vector_timeseries, FluorescenceTrace) and (angle_coords is None):
         if (
             isinstance(vector_timeseries.angle, np.ndarray) 
             and (vector_timeseries.angle.dtype == np.complex128)
@@ -156,7 +166,12 @@ def pva(
             and all(x is not None for x in vector_timeseries.angle)
         ):
             angle_coords = np.exp(-1j*vector_timeseries.angle)
+    elif angle_coords is None:
+        angle_coords = np.exp(np.linspace(np.pi, -np.pi, vector_timeseries.shape[0])*1j) # it goes clockwise.
 
+    if angle_coords.dtype != np.complex128:
+        angle_coords = np.exp(1j*angle_coords)
+
     pva_val = np.asarray(
         np.matmul(
             angle_coords,

siffpy/siffmath/phase/traces.py

@@ -172,7 +172,7 @@ def invert(self):
         """
         self[...] = 1.0/self
 
-    def inverted(self)->np.ndarray:
+    def inverted(self)->'PhaseTrace':
         """
         Returns a _new_ PhaseTrace object with the phase inverted without
         modifying the existing trace.