diff options
-rw-r--r-- | libcrystfel/src/peaks.c | 78 | ||||
-rw-r--r-- | libcrystfel/src/peaks.h | 3 | ||||
-rw-r--r-- | src/im-sandbox.c | 2 |
3 files changed, 83 insertions, 0 deletions
diff --git a/libcrystfel/src/peaks.c b/libcrystfel/src/peaks.c index 503d534d..bf6a92f0 100644 --- a/libcrystfel/src/peaks.c +++ b/libcrystfel/src/peaks.c @@ -831,3 +831,81 @@ void integrate_reflections(struct image *image, int use_closer, int bgsub, free(il); } + + +void validate_peaks(struct image *image, double min_snr, + int ir_inn, int ir_mid, int ir_out) +{ + int i, n; + ImageFeatureList *flist; + int n_wtf, n_int, n_dft, n_snr, n_prx; + + flist = image_feature_list_new(); + if ( flist == NULL ) return; + + n = image_feature_count(image->features); + + /* Loop over peaks, putting each one through the integrator */ + n_wtf = 0; n_int = 0; n_dft = 0; n_snr = 0; n_prx = 0; + for ( i=0; i<n; i++ ) { + + struct imagefeature *f; + int r; + double d; + int idx; + double f_fs, f_ss; + double intensity, sigma; + struct panel *p; + + f = image_get_feature(image->features, i); + if ( f == NULL ) { + n_wtf++; + continue; + } + + p = find_panel(image->det, f->fs, f->ss); + if ( p == NULL ) { + n_wtf++; + continue; + } + + r = integrate_peak(image, f->fs, f->ss, + &f_fs, &f_ss, &intensity, &sigma, + ir_inn, ir_mid, ir_out, 0, NULL, NULL); + if ( r ) { + n_int++; + continue; + } + + /* It is possible for the centroid to fall outside the image */ + if ( (f_fs < p->min_fs) || (f_fs > p->max_fs) + || (f_ss < p->min_ss) || (f_ss > p->max_ss) ) + { + n_dft++; + continue; + } + + if ( fabs(intensity)/sigma < min_snr ) { + n_snr++; + continue; + } + + /* Check for a nearby feature */ + image_feature_closest(flist, f_fs, f_ss, &d, &idx); + if ( d < 2.0*ir_inn ) { + n_prx++; + continue; + } + + /* Add using "better" coordinates */ + image_add_feature(flist, f_fs, f_ss, image, intensity, NULL); + + } + + //STATUS("HDF5: %i peaks, validated: %i. WTF: %i, integration: %i," + // " drifted: %i, SNR: %i, proximity: %i\n", + // n, image_feature_count(flist), + // n_wtf, n_int, n_dft, n_snr, n_prx); + image_feature_list_free(image->features); + image->features = flist; +} diff --git a/libcrystfel/src/peaks.h b/libcrystfel/src/peaks.h index afc7b868..39fdf54c 100644 --- a/libcrystfel/src/peaks.h +++ b/libcrystfel/src/peaks.h @@ -53,4 +53,7 @@ extern int peak_sanity_check(struct image *image); extern void estimate_resolution(RefList *list, UnitCell *cell, double *min, double *max); +extern void validate_peaks(struct image *image, double min_snr, + int ir_inn, int ir_mid, int ir_out); + #endif /* PEAKS_H */ diff --git a/src/im-sandbox.c b/src/im-sandbox.c index c313bda9..31dff291 100644 --- a/src/im-sandbox.c +++ b/src/im-sandbox.c @@ -291,6 +291,8 @@ static void process_image(const struct index_args *iargs, iargs->hdf5_peak_path)) { ERROR("Failed to get peaks from HDF5 file.\n"); } + validate_peaks(&image, iargs->min_int_snr, + iargs->ir_inn, iargs->ir_mid, iargs->ir_out); break; case PEAK_ZAEF: |