@@ -279,28 +279,68 @@ int getStartTimeInSet(const std::vector<eventTimeTrack>& tracks, std::vector<int
279279 return 2 ; // no event time in the set
280280 }
281281
282- unsigned long ncomb = combinatorial[ntracks];
283- for (unsigned long comb = 0 ; comb < ncomb; comb++) {
282+ // Pre-filter: for each track, determine which hypotheses pass the time cut
283+ double preStarttime[maxNtrackInSet][3 ];
284+ double preWeighttime[maxNtrackInSet][3 ];
285+ int validHypo[maxNtrackInSet][3 ]; // which hypothesis indices are valid
286+ int nValidHypo[maxNtrackInSet]; // how many valid hypotheses per track
287+
288+ unsigned long ncombReduced = 1 ;
289+ for (int itrk = 0 ; itrk < ntracks; itrk++) {
290+ const eventTimeTrack& ctrack = tracks[trackInSet[itrk]];
291+ nValidHypo[itrk] = 0 ;
292+ for (int h = 0 ; h < 3 ; h++) {
293+ double st = ctrack.mSignal - ctrack.expTimes [h];
294+ if (std::abs (st - refT0) < 2000 ) {
295+ int idx = nValidHypo[itrk];
296+ validHypo[itrk][idx] = h;
297+ preStarttime[itrk][idx] = st;
298+ preWeighttime[itrk][idx] = 1 . / (ctrack.expSigma [h] * ctrack.expSigma [h]);
299+ nValidHypo[itrk]++;
300+ }
301+ }
302+ if (nValidHypo[itrk] == 0 ) {
303+ // No valid hypothesis for this track; treat as 1 entry with zero weight
304+ // to keep combination enumeration simple
305+ validHypo[itrk][0 ] = 0 ;
306+ preStarttime[itrk][0 ] = 0 ;
307+ preWeighttime[itrk][0 ] = 0 ;
308+ nValidHypo[itrk] = 1 ;
309+ }
310+ ncombReduced *= nValidHypo[itrk];
311+ }
312+
313+ // Enumerate only valid hypothesis combinations
314+ for (unsigned long comb = 0 ; comb < ncombReduced; comb++) {
284315 unsigned long curr = comb;
285316
286317 int ngood = 0 ;
287318 double average = 0 ;
288319 double sumweights = 0 ;
289- // get track info in the set for current combination
320+ // Decode the combination using per-track valid hypothesis counts
321+ unsigned long origComb = 0 ;
322+ unsigned long origBase = 1 ;
290323 for (int itrk = 0 ; itrk < ntracks; itrk++) {
291- hypo[itrk] = curr % 3 ;
292- curr /= 3 ;
293- const eventTimeTrack& ctrack = tracks[trackInSet[itrk]];
294- starttime[itrk] = ctrack.mSignal - ctrack.expTimes [hypo[itrk]];
295-
296- if (std::abs (starttime[itrk] - refT0) < 2000 ) { // otherwise time inconsistent with the int BC
297- weighttime[itrk] = 1 . / (ctrack.expSigma [hypo[itrk]] * ctrack.expSigma [hypo[itrk]]);
324+ int localIdx = curr % nValidHypo[itrk];
325+ curr /= nValidHypo[itrk];
326+ int h = validHypo[itrk][localIdx];
327+ hypo[itrk] = h;
328+ origComb += h * origBase;
329+ origBase *= 3 ;
330+ starttime[itrk] = preStarttime[itrk][localIdx];
331+ weighttime[itrk] = preWeighttime[itrk][localIdx];
332+
333+ if (weighttime[itrk] > 0 ) {
298334 average += starttime[itrk] * weighttime[itrk];
299335 sumweights += weighttime[itrk];
300336 ngood++;
301337 }
302338 }
303339
340+ if (ngood < 2 ) {
341+ continue ;
342+ }
343+
304344 average /= sumweights;
305345
306346 // compute chi2
@@ -313,7 +353,7 @@ int getStartTimeInSet(const std::vector<eventTimeTrack>& tracks, std::vector<int
313353 chi2 /= (ngood - 1 );
314354
315355 if (chi2 < chi2best) {
316- bestComb = comb ;
356+ bestComb = origComb ;
317357 chi2best = chi2;
318358 averageBest = average;
319359 }
0 commit comments