diff options
Diffstat (limited to 'src/intensities.c')
-rw-r--r-- | src/intensities.c | 154 |
1 files changed, 73 insertions, 81 deletions
diff --git a/src/intensities.c b/src/intensities.c index f1bdf19..86f3de1 100644 --- a/src/intensities.c +++ b/src/intensities.c @@ -3,7 +3,7 @@ * * Extract integrated intensities by relrod estimation * - * (c) 2007-2008 Thomas White <taw27@cam.ac.uk> + * (c) 2007-2009 Thomas White <taw27@cam.ac.uk> * * dtr - Diffraction Tomography Reconstruction * @@ -24,8 +24,8 @@ /* Extract integrated reflection intensities by estimating the spike function * based on the observed intensity and the calculated excitation error from * the lattice refinement. Easy. */ -void intensities_extract(ControlContext *ctx) { - +void intensities_extract(ControlContext *ctx) +{ int i, j; int n_meas, n_dupl, n_notf; double max; @@ -36,123 +36,111 @@ void intensities_extract(ControlContext *ctx) { reflectionlist_free(ctx->integrated); } ctx->integrated = reflectionlist_new(); - + n_meas = 0; n_dupl = 0; n_notf = 0; max = 0; for ( i=0; i<ctx->images->n_images; i++ ) { - + ImageRecord *image; - + image = &ctx->images->images[i]; - if ( image->rflist == NULL ) image->rflist = reproject_get_reflections(image, ctx->cell_lattice, ctx); - + if ( image->rflist == NULL ) + image->rflist = reproject_get_reflections(image, + ctx->cell_lattice); + for ( j=0; j<image->rflist->n_features; j++ ) { - + ImageFeature *feature; signed int h, k, l; - + feature = &image->rflist->features[j]; - + h = feature->reflection->h; k = feature->reflection->k; l = feature->reflection->l; - + if ( feature->partner != NULL ) { - + if ( (h!=0) || (k!=0) || (l!=0) ) { - + double intensity; Reflection *ref; - - /* Perform relrod calculation of doom here. - * TODO: Figure out if this is even possible. */ + intensity = feature->partner->intensity; - - ref = reflectionlist_find(ctx->integrated, h, k, l); - + + ref = reflectionlist_find( + ctx->integrated, h, k, l); + if ( ref == NULL ) { - + Reflection *new; - - printf("IN: Adding %3i %3i %3i, intensity=%f\n", h, k, l, intensity); - - new = reflection_add(ctx->integrated, - feature->reflection->x, feature->reflection->y, feature->reflection->z, - intensity, REFLECTION_GENERATED); - + + new = reflection_add( + ctx->integrated, + feature->reflection->x, + feature->reflection->y, + feature->reflection->z, + intensity, + REFLECTION_GENERATED); + new->h = h; new->k = k; new->l = l; - - if ( intensity > max ) max = intensity; - + + if ( intensity > max ) + max = intensity; + n_meas++; - + } else { - - printf("IN: Duplicate measurement for %3i %3i %3i - ", h, k, l); - + if ( intensity > ref->intensity ) { - - printf("stronger.\n"); - + ref->x = feature->reflection->x; ref->y = feature->reflection->y; ref->z = feature->reflection->z; ref->intensity = intensity; - - } else { - printf("weaker.\n"); + } - + n_dupl++; - + } - + } - + } else { //printf("IN: %3i %3i %3i not found\n", h, k, l); n_notf++; } - + } - + } - + /* Normalise all reflections to the most intense reflection */ reflection = ctx->integrated->reflections; while ( reflection ) { - reflection->intensity /= max; - - /* Test mode */ - // if ( (reflection->h == 0) && (reflection->k == -1) && (reflection->l == 0) ) { - // reflection->intensity = 1.0; - // } else { - // reflection->intensity /= 10; - // } - reflection = reflection->next; - } - + printf("IN: %5i intensities measured\n", n_meas); printf("IN: %5i duplicated measurements\n", n_dupl); printf("IN: %5i predicted reflections not found\n", n_notf); - } -static int intensities_do_save(ReflectionList *integrated, Basis *cell, const char *filename) { - +static int intensities_do_save(ReflectionList *integrated, Basis *cell, + const char *filename) +{ FILE *fh; Reflection *reflection; UnitCell rcell; - + fh = fopen(filename, "w"); - + rcell = basis_get_cell(cell); fprintf(fh, "a %12.8f\n", rcell.a*1e9); fprintf(fh, "b %12.8f\n", rcell.b*1e9); @@ -160,46 +148,50 @@ static int intensities_do_save(ReflectionList *integrated, Basis *cell, const ch fprintf(fh, "alpha %12.8f\n", rad2deg(rcell.alpha)); fprintf(fh, "beta %12.8f\n", rad2deg(rcell.beta)); fprintf(fh, "gamma %12.8f\n", rad2deg(rcell.gamma)); - + reflection = integrated->reflections; while ( reflection ) { - fprintf(fh, "%3i %3i %3i %12.8f\n", reflection->h, reflection->k, reflection->l, reflection->intensity); + fprintf(fh, "%3i %3i %3i %12.8f\n", + reflection->h, reflection->k, reflection->l, + reflection->intensity); reflection = reflection->next; } fclose(fh); - - return 0; + return 0; } -static gint intensities_save_response(GtkWidget *widget, gint response, ControlContext *ctx) { - +static gint intensities_save_response(GtkWidget *widget, gint response, + ControlContext *ctx) +{ if ( response == GTK_RESPONSE_ACCEPT ) { char *filename; - filename = gtk_file_chooser_get_filename(GTK_FILE_CHOOSER(widget)); - if ( intensities_do_save(ctx->integrated, ctx->cell, filename) ) { - displaywindow_error("Failed to save cache file.", ctx->dw); + filename = gtk_file_chooser_get_filename( + GTK_FILE_CHOOSER(widget)); + if ( intensities_do_save(ctx->integrated, + ctx->cell, filename) ) { + displaywindow_error("Failed to save cache file.", + ctx->dw); } g_free(filename); } - + gtk_widget_destroy(widget); return 0; - } -void intensities_save(ControlContext *ctx) { - +void intensities_save(ControlContext *ctx) +{ GtkWidget *save; - save = gtk_file_chooser_dialog_new("Save Reflections to File", GTK_WINDOW(ctx->dw->window), + save = gtk_file_chooser_dialog_new("Save Reflections to File", + GTK_WINDOW(ctx->dw->window), GTK_FILE_CHOOSER_ACTION_SAVE, GTK_STOCK_CANCEL, GTK_RESPONSE_CANCEL, GTK_STOCK_SAVE, GTK_RESPONSE_ACCEPT, NULL); - g_signal_connect(G_OBJECT(save), "response", G_CALLBACK(intensities_save_response), ctx); + g_signal_connect(G_OBJECT(save), "response", + G_CALLBACK(intensities_save_response), ctx); gtk_widget_show_all(save); - } - |