aboutsummaryrefslogtreecommitdiff
path: root/src
diff options
context:
space:
mode:
Diffstat (limited to 'src')
-rw-r--r--src/rejection.c35
1 files changed, 26 insertions, 9 deletions
diff --git a/src/rejection.c b/src/rejection.c
index 7b78803e..9f298c0b 100644
--- a/src/rejection.c
+++ b/src/rejection.c
@@ -124,6 +124,7 @@ static double calculate_cchalf(RefList *template, RefList *full,
signed int h, k, l;
int nc = 0;
Reflection *refl;
+ double res;
get_indices(trefl, &h, &k, &l);
@@ -139,19 +140,29 @@ static double calculate_cchalf(RefList *template, RefList *full,
c = get_contributions(refl);
assert(c != NULL);
+ /* Calculate the resolution just once, using the cell from the
+ * first crystal to contribute, otherwise it takes too long */
+ res = resolution(crystal_get_cell(c->contrib_crystals[0]),
+ h, k, l);
+
/* Mean of contributions */
refl_sum = 0.0;
for ( j=0; j<c->n_contrib; j++ ) {
- double Ii;
+ double Ii, G, B;
if ( c->contrib_crystals[j] == exclude ) {
continue;
}
- /* FIXME: Apply corrections */
+ G = crystal_get_osf(c->contrib_crystals[j]);
+ B = crystal_get_Bfac(c->contrib_crystals[j]);
Ii = get_intensity(c->contribs[j]);
+ /* Total (multiplicative) correction factor */
+ Ii *= 1.0/G * exp(B*res*res) * get_lorentz(c->contribs[j])
+ / get_partiality(c->contribs[j]);
+
refl_sum += Ii;
nc++;
@@ -164,15 +175,20 @@ static double calculate_cchalf(RefList *template, RefList *full,
refl_sumsq = 0.0;
for ( j=0; j<c->n_contrib; j++ ) {
- double Ii;
+ double Ii, G, B;
if ( c->contrib_crystals[j] == exclude ) {
continue;
}
- /* FIXME: Apply corrections */
+ G = crystal_get_osf(c->contrib_crystals[j]);
+ B = crystal_get_Bfac(c->contrib_crystals[j]);
Ii = get_intensity(c->contribs[j]);
+ /* Total (multiplicative) correction factor */
+ Ii *= 1.0/G * exp(B*res*res) * get_lorentz(c->contribs[j])
+ / get_partiality(c->contribs[j]);
+
refl_sumsq += (Ii-refl_mean)*(Ii-refl_mean);
/* (nc already summed above) */
@@ -226,12 +242,13 @@ static void check_deltacchalf(Crystal **crystals, int n, RefList *full)
RefList *template = crystal_get_reflections(crystals[i]);
cchalf = calculate_cchalf(template, full, NULL, NULL);
cchalfi = calculate_cchalf(template, full, crystals[i], &nref);
- STATUS("Frame %i:", i);
- STATUS(" With = %f ", cchalf*100.0);
- STATUS("Without = %f", cchalfi*100.0);
- STATUS(" Delta = %f ", (cchalf - cchalfi)*100.0);
- STATUS("(nref = %i)\n", nref);
+ //STATUS("Frame %i:", i);
+ //STATUS(" With = %f ", cchalf*100.0);
+ //STATUS("Without = %f", cchalfi*100.0);
+ //STATUS(" Delta = %f ", (cchalf - cchalfi)*100.0);
+ //STATUS("(nref = %i)\n", nref);
vals[i] = cchalf - cchalfi;
+ progress_bar(i, n-1, "Calculating deltaCChalf");
}
mean = gsl_stats_mean(vals, 1, n);