diff options
| -rw-r--r-- | src/peaks.c | 10 | ||||
| -rw-r--r-- | src/process_hkl.c | 44 | ||||
| -rw-r--r-- | tests/integration_check.c | 1 |
3 files changed, 31 insertions, 24 deletions
diff --git a/src/peaks.c b/src/peaks.c index aa56da8b..7db2ccc9 100644 --- a/src/peaks.c +++ b/src/peaks.c @@ -257,16 +257,12 @@ int integrate_peak(struct image *image, int cfs, int css, if ( in_bad_region(image->det, *pfs, *pss) ) return 1; if ( sigma != NULL ) { - /* - * first term is standard deviation of background per pixel + /* First term is standard deviation of background per pixel * sqrt(pixel_counts) - increase of error for integrated value * sqrt(2) - increase of error for background subtraction */ - *sigma = sqrt( noise_meansq/noise_counts - (noise_mean * noise_mean)) - *sqrt(2.0*pixel_counts); - /* printf(" counts %d %d \n", noise_counts, pixel_counts); - printf(" intensity, bg, diff, %f, %f, %f \n", total, pixel_counts*noise_mean, total - pixel_counts*noise_mean); - printf(" sigma = %f \n", *sigma); */ + *sigma = sqrt(noise_meansq/noise_counts-(noise_mean*noise_mean)) + * sqrt(2.0*pixel_counts); } if ( pbg != NULL ) { diff --git a/src/process_hkl.c b/src/process_hkl.c index 10597337..4dfdaa22 100644 --- a/src/process_hkl.c +++ b/src/process_hkl.c @@ -4,7 +4,8 @@ * Assemble and process FEL Bragg intensities * * (c) 2006-2011 Thomas White <taw@physics.org> - * 2011 Andrew Martin + * (c) 2011 Andrew Martin <andrew.martin@desy.de> + * * Part of CrystFEL - crystallography with a FEL * */ @@ -171,12 +172,16 @@ static void merge_pattern(RefList *model, RefList *new, int max_only, } else if ( pass == 2 ) { double dev = get_sum_squared_dev(model_version); - double refl_esd = get_esd_intensity(refl); - - set_sum_squared_dev(model_version, - // dev + pow(refl_esd,2.0) ); - dev + pow(intensity, 2.0) ); - // dev + pow(intensity-model_int, 2.0)); + + /* Other ways of estimating the sigma are possible, + * choose from: + * dev += pow(get_esd_intensity(refl), 2.0); + * dev += pow(intensity, 2.0); + * But alter the other part of the calculation below + * as well. */ + dev += pow(intensity - model_int, 2.0); + + set_sum_squared_dev(model_version, dev); if ( hist_vals != NULL ) { int p = *hist_n; @@ -371,19 +376,24 @@ static void merge_all(FILE *fh, RefList *model, refl = next_refl(refl, iter) ) { double sum_squared_dev = get_sum_squared_dev(refl); - double intensity = get_intensity(refl); int red = get_redundancy(refl); int h, k, l; + double esd; get_indices(refl,&h,&k,&l); - - /* - set_esd_intensity(refl, - sqrt(sum_squared_dev)/(double)red); - */ - - set_esd_intensity(refl, - sqrt((sum_squared_dev/(double)red) - pow(intensity,2.0) )/(double)red); - + + /* Other ways of estimating the sigma are possible, + * such as: + * + * double intensity = get_intensity(refl); + * esd = sqrt( (sum_squared_dev/(double)red) + * - pow(intensity,2.0) ) ); + * + * But alter the other part of the calculation above + * as well. */ + esd = sqrt(sum_squared_dev)/(double)red; + + set_esd_intensity(refl, esd); + } } diff --git a/tests/integration_check.c b/tests/integration_check.c index 61124f09..82278fdc 100644 --- a/tests/integration_check.c +++ b/tests/integration_check.c @@ -4,6 +4,7 @@ * Check peak integration * * (c) 2011 Thomas White <taw@physics.org> + * (c) 2011 Andrew Martin <andrew.martin@desy.de> * * Part of CrystFEL - crystallography with a FEL * |