tracking::digitization::StripClusterer Class Reference

Nearest-neighbour clustering of fitted silicon-strip hits on a single sensor. More...

#include <StripClusterer.h>

Classes
struct	ClusterCandidate

Public Member Functions
	StripClusterer (double seed_threshold=4.0, double neighbor_threshold=3.0, double cluster_threshold=4.0, double noise_sigma_adc=5.0, double mean_time_ns=0.0, double time_window_ns=-1.0, double neighbor_delta_t_ns=-1.0, double max_chi2_ndf=-1.0)
std::vector< ClusterCandidate >	findClusters (const std::vector< ldmx::FittedSiStripHit > &hits) const
	Cluster a set of fitted strip hits from a single sensor layer.

Private Member Functions
bool	passesSeedCuts (const ldmx::FittedSiStripHit &h) const
bool	passesNeighborCuts (const ldmx::FittedSiStripHit &h, double cluster_weighted_t, double cluster_total_amp) const
double	hitNoise (const ldmx::FittedSiStripHit &h) const
	Per-strip noise RMS to use for h: the hit's own measured noise if set (> 0), otherwise the uniform noise passed to the constructor.

Private Attributes
double	seed_threshold_
double	neighbor_threshold_
double	cluster_threshold_
double	noise_sigma_adc_
double	mean_time_ns_
double	time_window_ns_
double	neighbor_delta_t_ns_
double	max_chi2_ndf_

Detailed Description

Nearest-neighbour clustering of fitted silicon-strip hits on a single sensor.

Ported from HPS NearestNeighborRMSClusterer (Java). Thresholds are all expressed in units of the per-strip noise RMS so they are independent of the absolute gain.

Algorithm

1. Build a channel→hit map. If two hits share a strip, keep the one with smaller |t0| (closest to the expected hit time).
2. Mark strips with amplitude/noise ≥ neighbor_threshold as clusterable.
3. Mark clusterable strips that also satisfy amplitude/noise ≥ seed_threshold and pass the timing and chi2 quality cuts as seeds.
4. For each seed (processed in strip-index order), grow a cluster by BFS:
   • Pop a strip from the "unchecked" queue and add it to the cluster.
   • Examine its two nearest neighbours (strip ± 1).
   • If a neighbour is still clusterable and its |t0 − cluster_weighted_t| is within neighbor_delta_t_ns, add it to the queue and mark it unavailable.
5. Accept the cluster if Σ(amp) / √Σ(noise²) > cluster_threshold.

Output

Each accepted cluster is reported as a ClusterCandidate: centroid_strip Charge-weighted mean strip index (fractional). total_amplitude Σ amplitude [ADC counts]. time_ns Amplitude-weighted mean hit time [ns]. sigma_strip Position uncertainty [strips]: charge-weighted RMS around the centroid (uses actual charge-sharing profile); floored at 1/√12 (single-strip binary limit) when RMS = 0. n_strips Number of strips in the cluster. layer_id Sensor layer (copied from the input hits). strip_ids Constituent strip indices (for diagnostics).

Convert to physical coordinates in the caller: local_u = centroid_strip * readout_pitch_mm sigma_u = sigma_strip * readout_pitch_mm

Definition at line 48 of file StripClusterer.h.

Constructor & Destructor Documentation

StripClusterer()

tracking::digitization::StripClusterer::StripClusterer	(	double	seed_threshold = 4.0,
		double	neighbor_threshold = 3.0,
		double	cluster_threshold = 4.0,
		double	noise_sigma_adc = 5.0,
		double	mean_time_ns = 0.0,
		double	time_window_ns = -1.0,
		double	neighbor_delta_t_ns = -1.0,
		double	max_chi2_ndf = -1.0 )

Parameters

seed_threshold	Minimum amplitude/noise to seed a cluster (default 4).
neighbor_threshold	Minimum amplitude/noise for a strip to join a cluster (default 3).
cluster_threshold	Minimum Σamp / √Σnoise² for the cluster to be kept (default 4).
noise_sigma_adc	Per-strip noise RMS [ADC counts].
mean_time_ns	Expected hit time for the seed timing cut [ns] (default 0).
time_window_ns	Half-width of seed timing window [ns]; ≤ 0 disables (default -1).
neighbor_delta_t_ns	Max |t0_neighbour − cluster_t| to join cluster [ns]; ≤ 0 disables (default -1).
max_chi2_ndf	Max chi2/ndf for a hit to be used; ≤ 0 disables (default -1).

Definition at line 11 of file StripClusterer.cxx.

    : seed_threshold_(seed_threshold),
      neighbor_threshold_(neighbor_threshold),
      cluster_threshold_(cluster_threshold),
      noise_sigma_adc_(noise_sigma_adc),
      mean_time_ns_(mean_time_ns),
      time_window_ns_(time_window_ns),
      neighbor_delta_t_ns_(neighbor_delta_t_ns),
      max_chi2_ndf_(max_chi2_ndf) {}

Member Function Documentation

findClusters()

std::vector< StripClusterer::ClusterCandidate > tracking::digitization::StripClusterer::findClusters

(

const std::vector< ldmx::FittedSiStripHit > &

hits

)

const

Cluster a set of fitted strip hits from a single sensor layer.

Parameters

hits	All FittedSiStripHits on one sensor (mixed layers are allowed but each hit's layer_id is preserved in the result).

Returns

List of cluster candidates, one per accepted cluster.

Definition at line 50 of file StripClusterer.cxx.

                                                       {
  // -------------------------------------------------------------------------
  // Build channel → hit map.
  // If two hits land on the same strip, keep the one with smaller |t0|.
  // -------------------------------------------------------------------------
  std::map<int, const ldmx::FittedSiStripHit*> channel_map;
  for (const auto& h : hits) {
    const int ch = h.getStripID();
    auto it = channel_map.find(ch);
    if (it == channel_map.end()) {
      channel_map[ch] = &h;
    } else {
      // Keep the hit closest to the expected hit time
      if (std::abs(h.getT0() - mean_time_ns_) <
          std::abs(it->second->getT0() - mean_time_ns_)) {
        it->second = &h;
      }
    }
  }
 
  // -------------------------------------------------------------------------
  // Determine which strips are clusterable (≥ neighbor threshold) and which
  // can seed a cluster (≥ seed threshold + timing/chi2 cuts).
  //
  // Thresholds are in units of the per-strip noise RMS.  Each hit may carry its
  // own measured noise (getNoise() > 0, from the real-data pedestal table);
  // when it does not (MC), we fall back to the uniform ctor noise so the MC
  // path is unchanged.
  // -------------------------------------------------------------------------
  std::set<int> clusterable_set;
  std::vector<int> seed_channels;
 
  for (const auto& [ch, hp] : channel_map) {
    const double amp = hp->getAmplitude();
    const double noise = hitNoise(*hp);
    if (amp >= neighbor_threshold_ * noise) {
      clusterable_set.insert(ch);
    }
    if (amp >= seed_threshold_ * noise && passesSeedCuts(*hp)) {
      seed_channels.push_back(ch);
    }
  }
 
  // Sort seeds by amplitude (highest first) so the strongest hit initiates.
  std::sort(seed_channels.begin(), seed_channels.end(), [&](int a, int b) {
    return channel_map.at(a)->getAmplitude() >
           channel_map.at(b)->getAmplitude();
  });
 
  // -------------------------------------------------------------------------
  // BFS expansion from each seed.
  // -------------------------------------------------------------------------
  std::vector<ClusterCandidate> clusters;
 
  for (int seed_ch : seed_channels) {
    // The seed might have already been claimed by an earlier cluster.
    if (clusterable_set.find(seed_ch) == clusterable_set.end()) continue;
 
    ClusterCandidate cand;
    double cluster_weighted_t = 0.0;
    double cluster_total_amp = 0.0;
    double cluster_noise_sq = 0.0;
 
    std::deque<int> unchecked;
    unchecked.push_back(seed_ch);
    clusterable_set.erase(seed_ch);
 
    while (!unchecked.empty()) {
      const int cur_ch = unchecked.front();
      unchecked.pop_front();
 
      const ldmx::FittedSiStripHit& hit = *channel_map.at(cur_ch);
      const double amp = hit.getAmplitude();
 
      // Accumulate cluster quantities.
      const double noise = hitNoise(hit);
      cand.strip_ids.push_back(cur_ch);
      cluster_total_amp += amp;
      cluster_weighted_t += amp * hit.getT0();
      cluster_noise_sq += noise * noise;
 
      // Check nearest neighbours (strip ± 1).
      for (int delta : {-1, +1}) {
        const int nb_ch = cur_ch + delta;
        if (clusterable_set.find(nb_ch) == clusterable_set.end()) continue;
 
        // Timing consistency with cluster so far.
        if (!passesNeighborCuts(*channel_map.at(nb_ch), cluster_weighted_t,
                                cluster_total_amp)) {
          continue;
        }
 
        unchecked.push_back(nb_ch);
        clusterable_set.erase(nb_ch);
      }
    }
 
    // Cluster S/N cut.
    if (cluster_noise_sq <= 0.0) continue;
    if (cluster_total_amp / std::sqrt(cluster_noise_sq) < cluster_threshold_) {
      continue;
    }
 
    // -----------------------------------------------------------------------
    // Compute the charge-weighted centroid strip and timing.
    // -----------------------------------------------------------------------
    double sum_amp_strip = 0.0;
    for (int ch : cand.strip_ids) {
      sum_amp_strip += channel_map.at(ch)->getAmplitude() * ch;
    }
 
    cand.centroid_strip = sum_amp_strip / cluster_total_amp;
    cand.total_amplitude = cluster_total_amp;
    cand.time_ns = cluster_weighted_t / cluster_total_amp;
    cand.n_strips = static_cast<int>(cand.strip_ids.size());
 
    // Charge-weighted RMS around the centroid [strips].
    // For multi-strip clusters this uses the actual charge-sharing profile,
    // giving sub-pitch resolution when charge is sharply peaked.
    // For single-strip clusters the RMS is zero, so we floor at the binary
    // single-strip uncertainty 1/√12.
    double sum_amp_dsq = 0.0;
    for (int ch : cand.strip_ids) {
      double d = ch - cand.centroid_strip;
      sum_amp_dsq += channel_map.at(ch)->getAmplitude() * d * d;
    }
    constexpr double k_single_strip_sigma = 1.0 / 3.4641;  // 1/√12
    const double rms = std::sqrt(sum_amp_dsq / cluster_total_amp);
    cand.sigma_strip = (rms > 0.0) ? rms : k_single_strip_sigma;
    cand.layer_id = channel_map.at(seed_ch)->getLayerID();
 
    clusters.push_back(std::move(cand));
  }
 
  return clusters;
}

References tracking::digitization::StripClusterer::ClusterCandidate::centroid_strip, ldmx::FittedSiStripHit::getAmplitude(), ldmx::FittedSiStripHit::getT0(), hitNoise(), tracking::digitization::StripClusterer::ClusterCandidate::sigma_strip, tracking::digitization::StripClusterer::ClusterCandidate::time_ns, and tracking::digitization::StripClusterer::ClusterCandidate::total_amplitude.

hitNoise()

double tracking::digitization::StripClusterer::hitNoise ( const ldmx::FittedSiStripHit & h ) const

inlineprivate

Per-strip noise RMS to use for h: the hit's own measured noise if set (> 0), otherwise the uniform noise passed to the constructor.

This is how real-data hits (with per-channel pedestal noise) and MC hits (uniform) share one clustering path.

Definition at line 106 of file StripClusterer.h.

                                                       {
    const double n = h.getNoise();
    return n > 0.0 ? n : noise_sigma_adc_;
  }

References ldmx::FittedSiStripHit::getNoise().

Referenced by findClusters().

passesNeighborCuts()

bool tracking::digitization::StripClusterer::passesNeighborCuts ( const ldmx::FittedSiStripHit & h, double cluster_weighted_t, double cluster_total_amp ) const

private

Definition at line 38 of file StripClusterer.cxx.

                                                                        {
  if (neighbor_delta_t_ns_ > 0.0 && cluster_total_amp > 0.0) {
    const double cluster_t = cluster_weighted_t / cluster_total_amp;
    if (std::abs(h.getT0() - cluster_t) > neighbor_delta_t_ns_) return false;
  }
  return true;
}

passesSeedCuts()

bool tracking::digitization::StripClusterer::passesSeedCuts ( const ldmx::FittedSiStripHit & h ) const

private

Definition at line 26 of file StripClusterer.cxx.

                                                                       {
  // Timing window (disabled if time_window_ns_ <= 0)
  if (time_window_ns_ > 0.0) {
    if (std::abs(h.getT0() - mean_time_ns_) > time_window_ns_) return false;
  }
  // Chi2/ndf quality cut (disabled if max_chi2_ndf_ <= 0)
  if (max_chi2_ndf_ > 0.0 && h.getNDF() > 0) {
    if (h.getReducedChi2() > max_chi2_ndf_) return false;
  }
  return true;
}

Member Data Documentation

cluster_threshold_

double tracking::digitization::StripClusterer::cluster_threshold_

private

Definition at line 113 of file StripClusterer.h.

max_chi2_ndf_

double tracking::digitization::StripClusterer::max_chi2_ndf_

private

Definition at line 118 of file StripClusterer.h.

mean_time_ns_

double tracking::digitization::StripClusterer::mean_time_ns_

private

Definition at line 115 of file StripClusterer.h.

neighbor_delta_t_ns_

double tracking::digitization::StripClusterer::neighbor_delta_t_ns_

private

Definition at line 117 of file StripClusterer.h.

neighbor_threshold_

double tracking::digitization::StripClusterer::neighbor_threshold_

private

Definition at line 112 of file StripClusterer.h.

noise_sigma_adc_

double tracking::digitization::StripClusterer::noise_sigma_adc_

private

Definition at line 114 of file StripClusterer.h.

seed_threshold_

double tracking::digitization::StripClusterer::seed_threshold_

private

Definition at line 111 of file StripClusterer.h.

time_window_ns_

double tracking::digitization::StripClusterer::time_window_ns_

private

Definition at line 116 of file StripClusterer.h.

---

The documentation for this class was generated from the following files:

• Tracking/include/Tracking/Digitization/StripClusterer.h
• Tracking/src/Tracking/Digitization/StripClusterer.cxx