[smeingast]

The group-level high-S/N error-floor correlations in scamp.xml
(AstromCorr_Internal_HighSN, AstromCorr_Reference_HighSN) are Pearson
coefficients but routinely come out > 1. In astrstats_fgroup()
(src/astrstats.c) they divide the full-sample numerator
(gcorr_int / gcorr_ref) by the high-S/N count, inflating the value by
~N_all / N_highSN; the correctly built high-S/N numerators
(gcorr_int_hsn / gcorr_ref_hsn, already scaled by the high-S/N sigma) were
never read — dead code. The fix uses them:

-                  gcorr_int_hsn += corr / (nmean_hsn - 1.0);
+                  gcorr_int_hsn += corr_hsn / (nmean_hsn - 1.0);
-  fgroup->sig_corr_int_hsn = nsource_hsn? gcorr_int / nsource_hsn : 0.0;
+  fgroup->sig_corr_int_hsn = nsource_hsn? gcorr_int_hsn / nsource_hsn : 0.0;
-  fgroup->sig_corr_ref_hsn = nmatch_hsn? gcorr_ref / nmatch_hsn : 0.0;
+  fgroup->sig_corr_ref_hsn = nmatch_hsn? gcorr_ref_hsn / nmatch_hsn : 0.0;

Diagnostic-only: no solver, cross-matching, or fitting code is touched.

Verification. A/B on a real 12-image VIRCAM stack, identical inputs, only
this patch differing: the two group values drop 2.92 → 0.53 and 2.02 → 0.26;
every other AstromCorr value (group non-HighSN and all per-field) is byte-for-byte
unchanged, and 0 values now exceed 1.

[smeingast]

Detailed root-cause analysis — kept out of the PR description for brevity. Expand below for the full write-up: the dead-code trace, the worked example with real FGroups counts, and why the non-HighSN group values and all per-field values are unaffected.

Full bug report (click to expand)

Group-level AstromCorr_*_HighSN in scamp.xml can exceed 1 (impossible for a Pearson correlation)

Summary

SCAMP writes group-level (FGroups) astrometric error-floor correlation
coefficients to scamp.xml:

• AstromCorr_Internal_HighSN
• AstromCorr_Reference_HighSN

These are Pearson correlation coefficients and must lie in [-1, 1]. In practice
the group-level high-S/N values routinely come out greater than 1. On a real
VIRCAM/VISIONS tile we observe a group AstromCorr_Reference_HighSN = 2.02 and
AstromCorr_Internal_HighSN = 2.92, while every per-field value and every
non-HighSN value in the same run is in range.

Affected version / commit

• Repository: https://github.com/astromatic/scamp
• Branch: master
• Commit at time of report: 2f608e1f598cea88880f7d8fcbfdc39e7a574e5e
• git describe: v2.15.0-1-g2f608e1
• configure.ac version string: 2.15.1dev
• The same logic is present in released 2.14.x binaries (the reference build
  used here reports SCAMP version 2.14.1).

File: src/astrstats.c, function astrstats_fgroup().

Symptom (real data)

From the FGroups table of an unpatched scamp.xml (12-image VIRCAM dither stack,
SN_THRESHOLDS 10.0,200.0):

AstromCorr_Internal               = 0.176884
AstromCorr_Internal_HighSN        = 2.92209     <-- impossible (> 1)
AstromCorr_Reference              = 0.117774
AstromCorr_Reference_HighSN       = 2.02049     <-- impossible (> 1)

The per-field (Fields table) values from the same run are all valid:

AstromCorr_Reference        (12 fields): 0.10 ... 0.27   (mean 0.187)
AstromCorr_Reference_HighSN (12 fields): 0.36 ... 0.60   (mean 0.4995)

So only the group-level *_HighSN numbers are broken.

Root cause

In astrstats_fgroup() the high-S/N group correlations are normalised with the
FULL-sample cross-product numerators (gcorr_int, gcorr_ref) but divided by the
HIGH-S/N counts (nsource_hsn, nmatch_hsn). The correctly built high-S/N
numerators (gcorr_int_hsn, gcorr_ref_hsn) are computed and even scaled by the
high-S/N sigma, then never read (dead code).

Line numbers below are from commit 2f608e1 (src/astrstats.c).

(a) Internal high-S/N numerator accumulated from the wrong cross-product (line 179)

Inside the high-S/N branch, corr_hsn is built about the high-S/N mean
(mean_hsn) at lines 173-178, but the accumulator adds the full-sample corr
instead:

160               corr = corr_hsn = 1.0;
161               for (i=0; i<naxis; i++)
162                 {
163                 dx = samp2->projpos[i] - mean[i];
164                 corr *= dx;
165                 lsig[i] += dx*dx / (nmean-1.0);
166                 }
167               gcorr_int += corr / (nmean - 1.0);
168               if (nmean_hsn>1)
169                 {
170                 sn2 = samp2->fluxerr>0.0? samp2->flux/samp2->fluxerr: 0.0;
171                 if (sn2 >= hsn_thresh)
172                   {
173                   for (i=0; i<naxis; i++)
174                     {
175                     dx = samp2->projpos[i] - mean_hsn[i];
176                     corr_hsn *= dx;
177                     lsig_hsn[i] += dx*dx / (nmean_hsn-1.0);
178                     }
179                   gcorr_int_hsn += corr / (nmean_hsn - 1.0);   // <-- uses corr, should be corr_hsn
180                   }
181                 }

(b) Internal high-S/N output uses the full-sample numerator (line 240)

After scaling by sigma (the high-S/N numerator is correctly scaled by the
high-S/N sigma at line 236, but that scaled value is never used):

234       gcorr_int *= fgroup->meanwcsscale[i] / fgroup->sig_interr[i];
...
236       gcorr_int_hsn *= fgroup->meanwcsscale[i] / fgroup->sig_interr_hsn[i];
...
239   fgroup->sig_corr_int = nsource? gcorr_int / nsource : 0.0;
240   fgroup->sig_corr_int_hsn = nsource_hsn? gcorr_int / nsource_hsn : 0.0;  // <-- gcorr_int, should be gcorr_int_hsn

(c) Reference high-S/N output uses the full-sample numerator (line 372)

Here the high-S/N numerator gcorr_ref_hsn is accumulated correctly (line 336
uses corr_hsn) and scaled by the high-S/N sigma (line 368), but the output line
still divides the full-sample gcorr_ref:

324       gcorr_ref += corr;                 // full sample
...
336         gcorr_ref_hsn += corr_hsn;       // high-S/N sample (correct)
...
366       gcorr_ref *= fgroup->meanwcsscale[i] / fgroup->sig_referr[i];
...
368       gcorr_ref_hsn *= fgroup->meanwcsscale[i] / fgroup->sig_referr_hsn[i];
...
371   fgroup->sig_corr_ref = nmatch? gcorr_ref / nmatch : 0.0;
372   fgroup->sig_corr_ref_hsn = nmatch_hsn? gcorr_ref / nmatch_hsn : 0.0;  // <-- gcorr_ref, should be gcorr_ref_hsn

Dead-code confirmation

gcorr_int_hsn and gcorr_ref_hsn are each declared, zeroed, accumulated, and
sigma-scaled, but neither is ever read into an fgroup-> output. Trace from the
source:

gcorr_int_hsn : declared 72, zeroed 109, accumulated 179, scaled 236  -> never read
gcorr_ref_hsn : declared 75, zeroed 249, accumulated 336, sub 362, scaled 368 -> never read

Why it inflates the value (and the worked example)

A Pearson correlation here is built as numerator / count, where the numerator
sums per-source cross-products each scaled by 1/sigma per axis. For the high-S/N
statistic the numerator must be summed over the high-S/N sample (N_highSN
terms) and divided by N_highSN; the buggy code substitutes the full-sample
numerator (a sum of N_all terms, scaled by the full-sample sigma) but still
divides by N_highSN. The numerator therefore contains roughly N_all positive
contributions while the denominator counts only N_highSN, inflating the result
by approximately N_all / N_highSN, modulated by the difference between the
full-sample and high-S/N sigma scalings.

Worked example from the real tile (FGroups counts in scamp.xml):

Reference:  AstromNDets_Reference = 27998 ,  AstromNDets_Reference_HighSN = 1632
            N_all / N_highSN = 17.2
Internal:   AstromNDets_Internal = 111665,  AstromNDets_Internal_HighSN = 6777
            N_all / N_highSN = 16.5

The full count ratio (~17x) over-states the visible inflation because the buggy
numerator is the full-sample one (scaled by the full-sample sigma, which is
larger than the high-S/N sigma, partially compensating), and corr/corr_hsn
are summed over different source sets. The net effect still pushes the value far
above 1:

                                 buggy (before)   fixed (after)
AstromCorr_Reference_HighSN          2.02049         0.258068
AstromCorr_Internal_HighSN           2.92209         0.530866

Both fixed values are now valid Pearson correlations (<= 1). The fixed group
reference value (0.26) is somewhat below the per-field high-S/N mean (0.50); this
is expected because the group statistic correlates the per-source mean residual
(averaged over the dither stack) against the reference, and averaging suppresses
correlated structure relative to the per-detection field statistic. The headline
result is that the impossible > 1 values disappear and the numbers become
physically interpretable.

Why field-level and non-HighSN values are unaffected

• The non-HighSN group values (sig_corr_int, sig_corr_ref) correctly pair the
  full-sample numerator (gcorr_int / gcorr_ref) with the full-sample count
  (nsource / nmatch), lines 239 and 371. They were always correct.
• All per-field correlations are computed in the separate per-field loop at the
  end of the function (lines 379-428), which uses each field's own accumulators
  (field->sig_corr_ref / field->sig_corr_ref_hsn and field->nchi2_ref*).
  That code never touches gcorr_int* / gcorr_ref*, so it is unaffected.

Fix

Three one-token changes in astrstats_fgroup():

1. Line 179: accumulate gcorr_int_hsn from corr_hsn (not corr).
2. Line 240: normalise sig_corr_int_hsn with gcorr_int_hsn (not gcorr_int).
3. Line 372: normalise sig_corr_ref_hsn with gcorr_ref_hsn (not gcorr_ref).

-                  gcorr_int_hsn += corr / (nmean_hsn - 1.0);
+                  gcorr_int_hsn += corr_hsn / (nmean_hsn - 1.0);

-  fgroup->sig_corr_int_hsn = nsource_hsn? gcorr_int / nsource_hsn : 0.0;
+  fgroup->sig_corr_int_hsn = nsource_hsn? gcorr_int_hsn / nsource_hsn : 0.0;

-  fgroup->sig_corr_ref_hsn = nmatch_hsn? gcorr_ref / nmatch_hsn : 0.0;
+  fgroup->sig_corr_ref_hsn = nmatch_hsn? gcorr_ref_hsn / nmatch_hsn : 0.0;

These make the high-S/N statistics use the high-S/N numerators (already scaled by
the high-S/N sigma) over the high-S/N counts, and they retire the previously dead
gcorr_int_hsn / gcorr_ref_hsn accumulators.

Verification

Built the patched tree (./autogen.sh && ./configure ... && make) on macOS arm64
and re-ran on the same 12 input LDAC catalogs and the same reference catalog and
config as the original pipeline invocation. Result: every AstromCorr_* value in
the patched scamp.xml satisfies |corr| <= 1 (0 violations), and the two
group *_HighSN values drop from 2.02 / 2.92 to 0.258 / 0.531. All non-HighSN
group values and all per-field values are byte-for-byte unchanged
(AstromCorr_Reference group 0.117774 before and after; per-field means 0.187
and 0.4995 before and after).

[ebertin]

Nice catch! Thank you @smeingast!