diff options
| -rw-r--r-- | src/hrs-scaling.c | 48 | ||||
| -rw-r--r-- | src/hrs-scaling.h | 5 | ||||
| -rw-r--r-- | src/partialator.c | 8 |
3 files changed, 54 insertions, 7 deletions
diff --git a/src/hrs-scaling.c b/src/hrs-scaling.c index 99490251..e2adaaf6 100644 --- a/src/hrs-scaling.c +++ b/src/hrs-scaling.c @@ -35,6 +35,45 @@ #define MAX_CYCLES (30) +char *find_common_reflections(struct image *images, int n) +{ + int i; + char *cref; + + cref = calloc(n*n, sizeof(char)); + + for ( i=0; i<n; i++ ) { + + Reflection *refl; + RefListIterator *iter; + + for ( refl = first_refl(images[i].reflections, &iter); + refl != NULL; + refl = next_refl(refl, iter) ) { + + signed int h, k, l; + int j; + + get_indices(refl, &h, &k, &l); + + for ( j=0; j<i; j++ ) { + Reflection *r2; + r2 = find_refl(images[j].reflections, h, k, l); + if ( r2 != NULL ) { + cref[i+n*j] = 1; + cref[j+n*i] = 1; + break; + } + } + + } + + } + + return cref; +} + + static void s_uhavha(signed int hat, signed int kat, signed int lat, struct image *image, double *uha, double *vha) { @@ -97,7 +136,7 @@ static double s_vh(struct image *images, int n, static double iterate_scale(struct image *images, int n, - ReflItemList *obs, const char *sym) + ReflItemList *obs, const char *sym, char *cref) { gsl_matrix *M; gsl_vector *v; @@ -181,6 +220,9 @@ static double iterate_scale(struct image *images, int n, /* Matrix is symmetric */ if ( b > a ) continue; + /* Zero if no common reflections */ + if ( cref[a+n*b] != 0 ) continue; + uhb = uha_arr[b]; vhb = vha_arr[b]; rhb = vhb - image_b->osf * uhb * Ih; @@ -336,7 +378,7 @@ static void normalise_osfs(struct image *images, int n) /* Scale the stack of images */ double *scale_intensities(struct image *images, int n, const char *sym, - ReflItemList *obs) + ReflItemList *obs, char *cref) { int m; double *I_full; @@ -350,7 +392,7 @@ double *scale_intensities(struct image *images, int n, const char *sym, i = 0; do { - max_shift = iterate_scale(images, n, obs, sym); + max_shift = iterate_scale(images, n, obs, sym, cref); STATUS("Iteration %2i: max shift = %5.2f\n", i, max_shift); i++; normalise_osfs(images, n); diff --git a/src/hrs-scaling.h b/src/hrs-scaling.h index ae384aa8..04a321e3 100644 --- a/src/hrs-scaling.h +++ b/src/hrs-scaling.h @@ -20,8 +20,9 @@ #include "image.h" -extern double *scale_intensities(struct image *image, int n, const char *sym, - ReflItemList *obs); +extern double *scale_intensities(struct image *images, int n, const char *sym, + ReflItemList *obs, char *cref); +extern char *find_common_reflections(struct image *images, int n); #endif /* HRS_SCALING_H */ diff --git a/src/partialator.c b/src/partialator.c index 309343cb..e99b2291 100644 --- a/src/partialator.c +++ b/src/partialator.c @@ -161,6 +161,7 @@ int main(int argc, char *argv[]) int n_found = 0; int n_expected = 0; int n_notfound = 0; + char *cref; /* Long options */ const struct option longopts[] = { @@ -379,12 +380,13 @@ int main(int argc, char *argv[]) STATUS("Mean measurements per unique reflection: %5.2f\n", (double)n_found / num_items(obs)); + cref = find_common_reflections(images, n_total_patterns); cts = new_list_count(); /* Make initial estimates */ STATUS("Performing initial scaling.\n"); select_scalable_reflections(images, n_total_patterns); - I_full = scale_intensities(images, n_total_patterns, sym, obs); + I_full = scale_intensities(images, n_total_patterns, sym, obs, cref); /* Iterate */ for ( i=0; i<n_iter; i++ ) { @@ -416,7 +418,8 @@ int main(int argc, char *argv[]) /* Re-estimate all the full intensities */ free(I_full); select_scalable_reflections(images, n_total_patterns); - I_full = scale_intensities(images, n_total_patterns, sym, obs); + I_full = scale_intensities(images, n_total_patterns, + sym, obs, cref); fclose(fhg); fclose(fhp); @@ -443,6 +446,7 @@ int main(int argc, char *argv[]) free(det); free(beam); free(cts); + free(cref); for ( i=0; i<n_total_patterns; i++ ) { cell_free(images[i].indexed_cell); free(images[i].filename); |