aboutsummaryrefslogtreecommitdiff
path: root/src/check_hkl.c
diff options
context:
space:
mode:
authorThomas White <taw@physics.org>2012-03-14 15:29:51 +0100
committerThomas White <taw@physics.org>2012-03-14 15:29:51 +0100
commit285cb2d55a9888083ab0224903f756c4c291c5a6 (patch)
tree0e3df97f5920dd0677ba900feca0f3c8c4fb236d /src/check_hkl.c
parenta52b5278dca146c4da556f4f93ec4863d26e459e (diff)
check_hkl: Fix inaccurate calculation of redundancy
Diffstat (limited to 'src/check_hkl.c')
-rw-r--r--src/check_hkl.c28
1 files changed, 21 insertions, 7 deletions
diff --git a/src/check_hkl.c b/src/check_hkl.c
index fa5fb0b7..649b86af 100644
--- a/src/check_hkl.c
+++ b/src/check_hkl.c
@@ -84,9 +84,10 @@ static void plot_shells(RefList *list, UnitCell *cell, const SymOpList *sym,
int i;
FILE *fh;
double snr_total = 0;
- int nmeas = 0;
+ int nrefl = 0;
int nmeastot = 0;
int nout = 0;
+ int nsilly = 0;
Reflection *refl;
RefListIterator *iter;
RefList *counted;
@@ -193,13 +194,15 @@ static void plot_shells(RefList *list, UnitCell *cell, const SymOpList *sym,
refl = next_refl(refl, iter) )
{
signed int h, k, l;
- double d;
+ double d, val, esd;
int bin;
int j;
get_indices(refl, &h, &k, &l);
d = resolution(cell, h, k, l) * 2.0;
+ val = get_intensity(refl);
+ esd = get_esd_intensity(refl);
bin = -1;
for ( j=0; j<NBINS; j++ ) {
@@ -210,6 +213,11 @@ static void plot_shells(RefList *list, UnitCell *cell, const SymOpList *sym,
}
if ( bin == -1 ) continue;
+ if ( !isfinite(val/esd) ) {
+ nsilly++;
+ continue;
+ }
+
measured[bin]++;
mean[bin] += get_intensity(refl);
}
@@ -253,20 +261,26 @@ static void plot_shells(RefList *list, UnitCell *cell, const SymOpList *sym,
measurements[bin] += get_redundancy(refl);
snr[bin] += val / esd;
snr_total += val / esd;
- nmeas++;
+ nrefl++;
nmeastot += get_redundancy(refl);
var[bin] += pow(val-mean[bin], 2.0);
}
- STATUS("overall <snr> = %f\n", snr_total/(double)nmeas);
+ STATUS("overall <snr> = %f\n", snr_total/(double)nrefl);
STATUS("%i measurements in total.\n", nmeastot);
- STATUS("%i reflections in total.\n", nmeas);
+ STATUS("%i reflections in total.\n", nrefl);
if ( nout ) {
STATUS("Warning; %i reflections outside resolution range.\n",
nout);
}
+ if ( nsilly ) {
+ STATUS("Warning; %i reflections had infinite or invalid values"
+ " of I/sigma(I).\n", nsilly);
+ }
+
+
fprintf(fh, "1/d centre # refs Possible Compl "
"Meas Red SNR Std dev Mean\n");
for ( i=0; i<NBINS; i++ ) {
@@ -278,9 +292,9 @@ static void plot_shells(RefList *list, UnitCell *cell, const SymOpList *sym,
cen*1.0e-9,
measured[i],
possible[i],
- 100.0*(float)measured[i]/possible[i],
+ 100.0*(double)measured[i]/possible[i],
measurements[i],
- (float)measurements[i]/measured[i],
+ (double)measurements[i]/measured[i],
snr[i]/(double)measured[i],
sqrt(var[i]/measured[i]),
mean[i]);