1. Phasor method for mono-exponential FLI Data processing#
A simulated, full 3-D fluorescence lifetime image dataset is generated as in the “Whole Image Simulation” example, and the phasor method is used to process the data.
This example walks through:
Repeating the whole-image simulation to generate the FLI data
The phasor computation
Phasor calibration with the IRF (if a known dye is used, it can also be used for calibration here — not shown in this example)
3-harmonics computation
Phasor plot (with pixel counts)
Mapping phasor color to pixel color (pixel-wise phasor overlay) and pixel-wise, intensity-weighted phasor color overlay
Phasor lifetime computation
Phasor plot of different harmonics
Curve fitting and fraction computation (based on user-provided lifetime values)
Author - Vikas
## importing the modules required for this work
import sys
import numpy as np
sys.path.insert(0, "ex_helper") # local helper for synthetic image generation
from sim_general_image import LettersShape, ROIMaskGenerator
from pyfli.analysis.utils import random_true_pixel
from pyfli.analyticalWorkflow import AnalyticalHelpers
from pyfli.data_cc import Normalization
from pyfli.data_text import MessageDisplay
from pyfli.data_vnp import ColorProcessor, DataViewer
from pyfli.io import DataOperations
from pyfli.phasor.phasorS import PhasorAnalyzer
from pyfli.simulator import FLIModelImageGenerator
1.1. Generating a test image#
A test image is simulated to compare the ground truth against the estimated values. In this example, a 128 x 512 image is generated with the letters “F”, “L”, “I”, “M”. Each letter is assigned a different intensity value (grayscale bit size) and a different lifetime: 0.7 ± 0.05 ns, 0.8 ± 0.05 ns, 0.9 ± 0.05 ns, and 1.0 ± 0.05 ns, respectively.
h, w = 128, 512
generator = ROIMaskGenerator((h, w), top=10, bottom=10, left=10, right=10)
FLIM_letters = LettersShape(letters=("F", "L", "I", "M"), gap=0.15)
custom_intensities = [0.7, 0.8, 0.9, 1.0]
# Preview the generated image
_ = generator.plot_preview(
FLIM_letters,
bit_depth=10,
intensities=custom_intensities,
show=True,
)
# Path to the local IRF file (e.g. an SPCImage-exported .txt calibration file)
IRF_PATH = "<select file/folder path>"
loader = DataOperations(irf_path=IRF_PATH)
irf_data = loader.load_irf()
gate_delay = 12.5 / irf_data.shape[2]
num_gates = irf_data.shape[2]
freq = AnalyticalHelpers(
laser_period=12.5, gate_delay=gate_delay, num_gate=num_gates
).freq_computation()
INFO:pyfli:Initiating IRF load from: <select file/folder path>
# select the type of decay model to simulate
MODEL_TYPE = "mono-exponential"
if MODEL_TYPE == "bi-exponential":
mono_fraction = 0.0
elif MODEL_TYPE == "mono-exponential":
mono_fraction = 1.0
else:
mono_fraction = 0.4 # fraction of mono-exponential model in sampled data
# Your constant, default configuration
BASE_CONFIG = {
# Modular Noise
"jitter": False, # offset artifact (due to jitter)
"dcr_on": True, # Dark Count Rate (thermal background)
"poisson": False, # Shot noise — Large detector only
"qe_on": True, # Quantum efficiency scaling
"read_noise_on": False, # Gaussian read noise — Large detector only
# Sensor
"sensor_type": "discrete", # "continuous" | "discrete"
"bit": 12,
"dcr": 0.08, # Mean dark counts per bin
"laser_feq": freq[1], # Laser repetition rate (MHz) → period = 12.5 ns
"round_on": True, # Round photon counts to integers — Macro_sim only
"clip_on": True, # Clip at bit-depth ceiling — Macro_sim: max_adc_val; TCSPC: max_bin_count
# Fluorescence Physics (FLI / FRET)
"tau2": (1, 1), # τ₂ ~ TruncNormal(mu=1 ns, sigma=0.5 ns)
"tau2_dist": "beta", # if dist "beta" or "normal" (default)
"efficiency": (1, 1), # FRET efficiency E ~ Beta(2, 5) → [0.1, 1.0]
"A1_fraction": (1, 1), # Amplitude fraction A₁ ~ Beta(2, 5) → [0.05, 0.95]
"photo_count": (2, 5), # Peak intensity ~ Beta(2, 5) × max_adc - large detectors
"mono_fraction": mono_fraction, # Fraction of pixels forced mono-exponential (0.0 = all bi-exp)
"n_cycles": (1_500_000, 2_000_000), # Accumulation cycles — for photon counter
}
Different ROIs (letter areas, in this case) are assigned different lifetimes.
def get_config(**kwargs):
"""Create a config by overriding defaults with specific test parameters."""
config = BASE_CONFIG.copy()
config.update(kwargs)
return config
ROI0 = {}
ROI1 = get_config(
tau2_beta_range=(0.05, 0.5),
)
ROI2 = get_config(
tau2_beta_range=(0.05, 0.7),
)
ROI3 = get_config(
tau2_beta_range=(0.05, 0.9),
)
ROI4 = get_config(
tau2_beta_range=(0.05, 1.1),
)
img_intensity = generator.generate_intensity_image(
FLIM_letters, intensities=custom_intensities
)
img_cluster = generator.generate_cluster_mask(FLIM_letters)
b_bool_mask = generator.generate_binary_mask(FLIM_letters)
simulated_img = FLIModelImageGenerator(
irf_data=irf_data[120, 40, :],
intensity_image=img_intensity,
roi_mask=img_cluster,
roi_params=[ROI0, ROI1, ROI2, ROI3, ROI4],
method="PHOTON_COUNTER",
verbose=True,
bool_mask=b_bool_mask,
)
gt_data = simulated_img.generate_image()
INFO:pyfli:Generating PHOTON_COUNTER FLI Image [128x512x256]...
res = gt_data["results"]["maps"]
if MODEL_TYPE == "bi-exponential":
data_list = [
res["tau1_map"],
res["tau2_map"],
res["alpha1_map"],
res["A1_map"],
res["A2_map"],
res["fret_efficiency_map"],
res["tau_mean_map"],
res["photon_count_map"],
res["mono_map"],
]
data_names = [
"tau1_map",
"tau2_map",
"alpha1_map",
"A1_map",
"A2_map",
"fret_efficiency_map",
"tau_mean_map",
"photon_count_map",
"mono_map",
]
rows = 3
fig_size = (13, 9)
else:
data_list = [res["tau_map"], res["photon_count_map"]]
data_names = ["tau_map", "photon_count_map"]
rows = 1
fig_size = (12, 3)
jet_m = ColorProcessor().lowest_zero("jet")
cmaps = [jet_m] * len(data_list)
v_ranges = None
px = None
if px is None:
cols = int(len(data_list) / rows)
else:
cols = int(len(data_list) / rows) + 1
_ = DataViewer().display_data(
data_list,
structure=(rows, cols),
coord=px,
data_names=data_names,
cmaps=cmaps,
v_ranges=v_ranges,
figsize=fig_size,
normalize=False,
yscale="linear",
)
Check the decay, IRF, and fit at a randomly selected pixel.
_decay = gt_data["raw_data"]["decay"]
_irf = gt_data["raw_data"]["irf"]
x, y = random_true_pixel(b_bool_mask)
TRs = gt_data["results"]["TR_maps"]
maps = gt_data["results"]["maps"]
print(f"the non-zero pixel selected for probing is ({x}, {y})")
irf_norm = Normalization(_irf).norm_scale(_decay)
DataViewer().plot_fli_px(
data_list=[_decay, irf_norm, TRs["fit_map"], TRs["residual_map"]],
pixel=(x, y),
mode=[0, 1, 2],
mode2=[1],
names=["decay", "irf", "fit"],
cmap=jet_m,
)
_ = MessageDisplay().get_pixel_summary(data_maps=maps, px=(x, y))
the non-zero pixel selected for probing is (44, 324)
INFO:pyfli:
Pixel (44, 324)
─────────────────────
A 11428.3057
α —
τ₁ 0.9241
τ₂ —
R² —
Red.χ² —
Raw.χ² —
v-shift —
h-shift —
─────────────────────
1.1.1. Phasor computation and visualization#
The framework below computes the phasor, harmonics, and phasor-based decay generation.
colorset = "rainbow"
freq_hz = freq[1] * 1e6
time_axis_ns = np.linspace(0, gate_delay * num_gates, _irf.shape[2])
phasor = PhasorAnalyzer(frequency_hz=freq_hz, time_axis_ns=time_axis_ns, n_harmonics=3)
G, S = phasor.create_phasor_gpu(_decay)
Gc, Sc = phasor.calibrate_pixelwise(G, S, _irf) # pixel-wise IRF calibration
tau_map_ns = phasor.compute_lifetime(Gc[0], Sc[0])
im1 = phasor.plot_phasor_diagram(
Gc[0],
Sc[0],
mask=b_bool_mask,
# colors="jet",
hexbin_color=colorset,
figsize=(7, 3),
half_circle=True,
kdeplot=False,
kde_color="red",
kde_levels=1,
kde_linewidths=1,
kde_alpha=0.5,
)
bg_color = "black" # either "black" or "white"
ph_cols_ = "viridis_r"
im2 = phasor.plot_overlay_subplots(
_decay,
Gc[0],
Sc[0],
mask=b_bool_mask,
colormaps=["plasma", "jet", "turbo_r", ph_cols_],
figsize=(13, 7),
xlim=(0, 1.1),
bg_color=bg_color,
transpose=False,
kdeplot=False,
)
1.2. Additional steps to plot different harmonics#
im3 = phasor.plot_phasor_harmonics(
Gc,
Sc,
harmonics=(1, 2, 3),
mask=b_bool_mask,
hexbin_color=colorset,
figsize=(15, 3),
)
1.2.1. Phasor of fitting based on the guessed lifetime#
Selecting the mono-exponential model for reconstructing the fit generates a general set of fit values.
x, y = random_true_pixel(b_bool_mask) # e.g. (62, 47)
fitting_type = "mono-exponential" # "mono-exponential" or "bi-exponential"
if fitting_type == "mono-exponential":
# user-provided input: tau_ns and pixel (x, y)
tau_ns = 0.9
phasor.plot_pixel_fit_single_exp(_irf, _decay, tau_map_ns, y, x, log_scale=False)
else:
tau1_ns, tau2_ns = 0.8, 1.3
A1, A2 = phasor.compute_fractions(Gc[0], Sc[0], tau1_ns, tau2_ns)
reconstructed_decay = phasor.analyze_biexponential_and_reconstruct(
Gc[0], Sc[0], _irf, tau1_ns=tau1_ns, tau2_ns=tau2_ns, plot=True
)
print(f"Reconstructed decay shape: {reconstructed_decay.shape}")
# Check whether the selected lifetime values make sense
phasor.plot_pixel_fit(_irf, _decay, reconstructed_decay, y, x, log_scale=False)
Selecting the bi-exponential model for reconstructing the fit generates a general set of fit values.
x, y = random_true_pixel(b_bool_mask)
fitting_type = "bi-exponential" # "mono-exponential" or "bi-exponential"
if fitting_type == "mono-exponential":
# user-provided input: tau_ns and pixel (x, y)
tau_ns = 0.9
phasor.plot_pixel_fit_single_exp(_irf, _decay, tau_map_ns, y, x, log_scale=False)
else:
tau1_ns, tau2_ns = 0.8, 1.3
A1, A2 = phasor.compute_fractions(Gc[0], Sc[0], tau1_ns, tau2_ns)
reconstructed_decay = phasor.analyze_biexponential_and_reconstruct(
Gc[0], Sc[0], _irf, tau1_ns=tau1_ns, tau2_ns=tau2_ns, plot=True
)
print(f"Reconstructed decay shape: {reconstructed_decay.shape}")
# Check whether the selected lifetime values make sense
phasor.plot_pixel_fit(_irf, _decay, reconstructed_decay, y, x, log_scale=False)
Reconstructed decay shape: (128, 512, 256)
