aboutsummaryrefslogtreecommitdiff
path: root/src
diff options
context:
space:
mode:
authorThomas White <taw@physics.org>2013-07-25 12:04:01 +0200
committerThomas White <taw@physics.org>2013-07-25 12:04:01 +0200
commit0431728cc5efd3321b1c60f97a830dd525cf04c8 (patch)
tree028cfdde49c4089fb09ce187a59fe74d98ef8bcf /src
parent0bb8df9b6a54db5da84d9fc9c26efcf5974fc89b (diff)
partialator: Simplify scaling
Diffstat (limited to 'src')
-rw-r--r--src/hrs-scaling.c38
-rw-r--r--src/partialator.c3
2 files changed, 23 insertions, 18 deletions
diff --git a/src/hrs-scaling.c b/src/hrs-scaling.c
index c249ef46..47c4301e 100644
--- a/src/hrs-scaling.c
+++ b/src/hrs-scaling.c
@@ -48,10 +48,7 @@
/* Maximum number of iterations of scaling per macrocycle. */
-#define MAX_CYCLES (50)
-
-/* ESD of restraint driving scale factors to unity */
-#define SCALING_RESTRAINT (1.0)
+#define MAX_CYCLES (10)
struct scale_queue_args
@@ -110,9 +107,8 @@ static void run_scale_job(void *vwargs, int cookie)
refl = next_refl(refl, iter) )
{
signed int h, k, l;
- double Ih, Ihl, esd;
+ double Ih, Ihl, esd, corr;
Reflection *r;
- double corr;
if ( !get_scalable(refl) ) continue;
@@ -151,19 +147,18 @@ static void run_scale_job(void *vwargs, int cookie)
Ihl = get_intensity(refl) / corr;
esd = get_esd_intensity(refl) / corr;
- num += Ih * (Ihl/G) / pow(esd/G, 2.0);
- den += pow(Ih, 2.0)/pow(esd/G, 2.0);
+ num += Ih * Ihl;
+ den += Ihl * Ihl;
}
- //num += image->osf / pow(SCALING_RESTRAINT, 2.0);
- //den += pow(image->osf, 2.0)/pow(SCALING_RESTRAINT, 2.0);
-
g = num / den;
if ( !isnan(g) && !isinf(g) ) {
- crystal_set_osf(cr, g*G);
+ crystal_set_osf(cr, g);
+ wargs->shift = fabs((g/G)-1.0);
+ } else {
+ wargs->shift = 0.0;
}
- wargs->shift = fabs(g-1.0);
}
@@ -298,8 +293,8 @@ static void run_merge_job(void *vwargs, int cookie)
Ihl = get_intensity(refl) / corr;
esd = get_esd_intensity(refl) / corr;
- num += (Ihl/G) / pow(esd/G, 2.0);
- den += 1.0 / pow(esd/G, 2.0);
+ num += Ihl * G;
+ den += 1.0;
red++;
set_temp1(f, num);
@@ -518,6 +513,7 @@ RefList *scale_intensities(Crystal **crystals, int n, RefList *gref,
i = 0;
do {
+ int j;
RefList *reference;
/* Refine against reference or current "full" estimates */
@@ -529,8 +525,12 @@ RefList *scale_intensities(Crystal **crystals, int n, RefList *gref,
max_corr = iterate_scale(crystals, n, reference, n_threads,
pmodel);
- //STATUS("Scaling iteration %2i: max correction = %5.2f\n",
- // i+1, max_corr);
+ STATUS("Scaling iteration %2i: max correction = %5.2f\n",
+ i+1, max_corr);
+ for ( j=0; j<10; j++ ) {
+ printf(" %5.2f", crystal_get_osf(crystals[j]));
+ }
+ printf("\n");
/* No reference -> generate list for next iteration */
if ( gref == NULL ) {
@@ -542,6 +542,10 @@ RefList *scale_intensities(Crystal **crystals, int n, RefList *gref,
} while ( (max_corr > 0.01) && (i < MAX_CYCLES) );
+ if ( i == MAX_CYCLES ) {
+ ERROR("Warning: Scaling did not converge.\n");
+ }
+
if ( gref != NULL ) {
full = lsq_intensities(crystals, n, n_threads, pmodel);
} /* else we already did it */
diff --git a/src/partialator.c b/src/partialator.c
index 73046e0a..07d7437c 100644
--- a/src/partialator.c
+++ b/src/partialator.c
@@ -190,7 +190,7 @@ static int select_scalable_reflections(RefList *list, RefList *reference)
if ( get_partiality(refl) < 0.1 ) sc = 0;
v = fabs(get_intensity(refl));
esd = get_esd_intensity(refl);
- if ( v < 0.5*esd ) sc = 0;
+ //if ( v < 0.5*esd ) sc = 0;
/* If we are scaling against a reference set, we additionally
* require that this reflection is in the reference list. */
@@ -497,6 +497,7 @@ int main(int argc, char *argv[])
n_crystals = 0;
images = NULL;
crystals = NULL;
+
do {
RefList *as;