Add tail-fit model and fix IRF/fit-window gaps from user testing

Tail fit: bare exponentials past the peak, ignoring the IRF, matching
the LAS X n-exponential tail fit. Selectable as --fit-model tail with
an optional free t0.

Fixes from NUS testing:
- validate_installation.py failed Tile Stitching and Complete Workflow
  because mock_data imported ptu_writer as a top-level module, which
  only resolves under pytest. Now 10/10 when run as a script.
- .pck IRF files were dispatched only from fit_cli, so the GUI and
  interactive paths sent them to the scatter-PTU reader. Added
  irf_from_measured_file and routed all three entry points through it.
- The IRF browse dialog filtered on FLIM data extensions plus xlsx, so
  .pck files were never listed.
- Fit window and exclusion bands were command-line only. Added start,
  end and exclusion fields to Expert Settings, validated on confirm,
  wired into the FOV and stitch paths.

Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>

Author: alex1075

fit_cli.py

@@ -5,8 +5,8 @@
 import matplotlib
 import argparse
 from flimkit.formats import FLIMFile
-from flimkit.FLIM.irf_tools import gaussian_irf_from_fwhm, irf_from_scatter_ptu, irf_from_pck, irf_from_xlsx, irf_from_xlsx_analytical, estimate_irf_from_decay_parametric, estimate_irf_from_decay_raw, reconstruct_irf_from_decay, compare_irfs, machine_irf_prompt
-from flimkit.FLIM.fitters import fit_summed, fit_per_pixel, MIN_PHOTONS_PERPIX
+from flimkit.FLIM.irf_tools import gaussian_irf_from_fwhm, irf_from_scatter_ptu, irf_from_pck, irf_from_measured_file, irf_from_xlsx, irf_from_xlsx_analytical, estimate_irf_from_decay_parametric, estimate_irf_from_decay_raw, reconstruct_irf_from_decay, compare_irfs, machine_irf_prompt
+from flimkit.FLIM.fitters import fit_summed, fit_per_pixel, fit_summed_tail, MIN_PHOTONS_PERPIX
 from flimkit.utils.plotting import plot_summed, plot_pixel_maps, plot_lifetime_histogram
 from flimkit.utils.misc import print_summary
 from flimkit.utils.xlsx_tools import load_xlsx
@@ -52,6 +52,14 @@ def single_FOV_flim_fit_cli():
                     help='IRF FWHM in ns. Default: 1 bin width from PTU '
                          '(e.g. 0.097 ns for 97 ps bins). Override for other systems.')
     ap.add_argument('--nexp',     type=int,   default=n_exp, choices=[1, 2, 3])
+    ap.add_argument('--fit-model', choices=['reconv', 'tail'], default='reconv',
+                    help="'reconv' (default) convolves the model with the IRF. "
+                         "'tail' fits bare exponentials past the peak and ignores "
+                         'the IRF entirely (Leica LAS X n-exponential tail fit).')
+    ap.add_argument('--fit-t0', action='store_true',
+                    help='Tail fit only: let the lifetime offset t0 float instead of '
+                         'pinning it to the decay peak. Correlated with the amplitudes, '
+                         'so leave it off unless the amplitudes matter more than stability.')
     ap.add_argument('--tau-min',  type=float, default=Tau_min, help='ns')
     ap.add_argument('--tau-max',  type=float, default=Tau_max, help='ns')
     ap.add_argument('--mode',     default=D_mode,
@@ -180,10 +188,7 @@ def single_FOV_flim_fit_cli():
     print(f"\n[4] Building IRF")
     sigma_max = MACHINE_IRF_SIGMA_MAX_FULL
     if args.irf is not None:
-        if str(args.irf).lower().endswith('.pck'):
-            irf_prompt = irf_from_pck(args.irf, ptu.n_bins, channel=args.channel)
-        else:
-            irf_prompt = irf_from_scatter_ptu(args.irf, ptu, channel=args.channel)
+        irf_prompt = irf_from_measured_file(args.irf, ptu, channel=args.channel)
         strategy = 'scatter_ptu'
         has_tail = False
         fit_sigma = False
@@ -268,7 +273,25 @@ def single_FOV_flim_fit_cli():
     fit_tvb = tvb_profile is not None
     global_popt = None
     global_summary = None
+    _tail_fit = args.fit_model == 'tail'
     def _run_summed():
+        if _tail_fit:
+            return fit_summed_tail(
+                decay, ptu.tcspc_res, ptu.n_bins,
+                fit_bg, args.nexp, args.tau_min, args.tau_max,
+                fit_t0=args.fit_t0,
+                optimizer=args.optimizer,
+                n_restarts=args.restarts,
+                de_popsize=args.de_population,
+                de_maxiter=args.de_maxiter,
+                workers=args.workers,
+                polish=not args.no_polish,
+                cost_function=args.cost_function,
+                tvb_profile=tvb_profile, fit_tvb=fit_tvb,
+                fit_start_ns=args.fit_start_ns,
+                fit_end_ns=args.fit_end_ns,
+                exclude_ns=parse_exclude_ns(args.exclude_ns),
+            )
         return fit_summed(
             decay, ptu.tcspc_res, ptu.n_bins,
             irf_prompt, has_tail, fit_bg, fit_sigma,
@@ -287,7 +310,8 @@ def _run_summed():
             exclude_ns=parse_exclude_ns(args.exclude_ns),
         )
     if args.mode in ('summed', 'both'):
-        print(f"\n[5] Summed decay fit  ({args.nexp}-exp, optimizer={args.optimizer})")
+        print(f"\n[5] Summed decay fit  ({args.nexp}-exp"
+              f"{' tail' if _tail_fit else ''}, optimizer={args.optimizer})")
         global_popt, global_summary = _run_summed()
         print_summary(global_summary, strategy, args.nexp)
         if not args.no_plots:
@@ -327,6 +351,9 @@ def _run_summed():
             fit_idx=global_summary.get('fit_idx'),
             free_tau=getattr(args, 'free_tau', False),
             tvb_profile=tvb_profile, fit_tvb=fit_tvb,
+            fit_model=args.fit_model,
+            fit_t0=args.fit_t0 if _tail_fit else False,
+            t0_fixed=(global_summary['t0_ns'] * 1e-9) if _tail_fit else 0.0,
         )
         roi_name = Path(args.ptu).stem
         save_weighted_tau_images(

flimkit/FLIM/fit_tools.py

@@ -70,6 +70,51 @@ def find_fit_start(decay: np.ndarray, irf_prompt: np.ndarray,
     return fit_start
 
 
+def find_tail_fit_start(decay, peak_bin, n_bins, frac=0.9):
+    d = np.asarray(decay, dtype=float)
+    peak_val = float(d[peak_bin]) if 0 <= peak_bin < n_bins else 0.0
+    if peak_val <= 0:
+        return min(peak_bin + 1, n_bins - 1)
+    below = np.where(d[peak_bin:] < frac * peak_val)[0]
+    if len(below) == 0:
+        return min(peak_bin + 1, n_bins - 1)
+    return int(min(peak_bin + below[0], n_bins - 1))
+
+def _build_bounds_tail(n_exp, tau_min, tau_max, decay_peak, fit_bg,
+                       bg_init=0.0, bg_upper=None,
+                       fit_t0=False, t0_init=0.0, t0_range=0.0,
+                       fit_tvb=False, tvb_init=0.0, tvb_upper=None):
+    lo = [tau_min] * n_exp + [0.0] * n_exp
+    hi = [tau_max] * n_exp + [10 * decay_peak] * n_exp
+    if fit_t0:
+        lo += [t0_init - t0_range];  hi += [t0_init + t0_range]
+    if fit_bg:
+        _bg_hi = bg_upper if bg_upper is not None else bg_init * 1.5 + 10
+        lo += [0.0];  hi += [_bg_hi]
+    if fit_tvb:
+        _tvb_hi = tvb_upper if tvb_upper is not None else tvb_init * 1.5 + 10
+        lo += [0.0];  hi += [_tvb_hi]
+    return lo, hi
+
+def _pack_p0_tail(n_exp, tau_min, tau_max, decay_peak, fit_bg, bg_init,
+                  tau_override=None, fit_t0=False, t0_init=0.0,
+                  fit_tvb=False, tvb_init=0.0):
+    if tau_override is not None:
+        taus0 = np.asarray(tau_override)
+    else:
+        tmin = max(tau_min, 1e-14) * 1.001
+        tmax = tau_max * 0.999
+        taus0 = np.logspace(np.log10(tmin), np.log10(tmax), n_exp)
+    amps0 = np.full(n_exp, decay_peak / n_exp)
+    base = np.concatenate([taus0, amps0])
+    if fit_t0:
+        base = np.concatenate([base, [t0_init]])
+    if fit_bg:
+        base = np.concatenate([base, [bg_init]])
+    if fit_tvb:
+        base = np.concatenate([base, [tvb_init]])
+    return base
+
 def find_fit_end(decay, peak_bin, tau_max_s, tcspc_res, n_bins):
     candidate = min(n_bins, peak_bin + int(6.0 * tau_max_s / tcspc_res))
     search_start = int(0.82 * n_bins)