From d0cf9ad6d2025cc3ac59e58511dde116d999b7e3 Mon Sep 17 00:00:00 2001 From: taw27 Date: Tue, 22 Jan 2008 12:43:28 +0000 Subject: Nooooo! Use RADIANS internally... git-svn-id: svn://cook.msm.cam.ac.uk:745/diff-tomo/dtr@260 bf6ca9ba-c028-0410-8290-897cf20841d1 --- src/displaywindow.c | 3 ++- src/glbits.c | 3 ++- src/image.h | 6 +++--- src/imagedisplay.c | 4 ++-- src/iprtest.c | 37 ++++++++----------------------------- src/mapping.c | 8 +++++--- src/mrc.c | 4 ++-- src/qdrp.c | 11 ++++++----- src/refine.c | 5 +++++ src/reproject.c | 4 ++-- 10 files changed, 37 insertions(+), 48 deletions(-) diff --git a/src/displaywindow.c b/src/displaywindow.c index 434508c..af17934 100644 --- a/src/displaywindow.c +++ b/src/displaywindow.c @@ -38,6 +38,7 @@ #include "imagedisplay.h" #include "intensities.h" #include "glbits.h" +#include "utils.h" static gint displaywindow_changeview(GtkWidget *widget, GtkRadioAction *action, DisplayWindow *dw) { @@ -242,7 +243,7 @@ static gint displaywindow_setaxis_response(GtkWidget *widget, gint response, Dis offset = gtk_entry_get_text(GTK_ENTRY(dw->tiltaxis_entry)); sscanf(offset, "%f", &off); - mapping_adjust_axis(dw->ctx, off); + mapping_adjust_axis(dw->ctx, deg2rad(off)); } diff --git a/src/glbits.c b/src/glbits.c index bf523fa..5aed16b 100644 --- a/src/glbits.c +++ b/src/glbits.c @@ -27,6 +27,7 @@ #include "trackball.h" #include "reflections.h" #include "image.h" +#include "utils.h" /* Utility function to load and compile a shader, checking the info log */ static GLhandleARB glbits_load_shader(const char *filename, GLenum type) { @@ -588,7 +589,7 @@ void glbits_prepare(DisplayWindow *dw) { * so rotate tilt axis anticlockwise by omega. * Since the rotation is about +z, this is already anticlockwise * when looking down z. */ - glRotatef(ctx->images->images[0].omega, 0.0, 0.0, 1.0); + glRotatef(rad2deg(ctx->images->images[0].omega), 0.0, 0.0, 1.0); glScalef(50.0, 1.0, 1.0); glTranslatef(-0.5, 0.0, 0.0); glMaterialfv(GL_FRONT_AND_BACK, GL_EMISSION, yellow); diff --git a/src/image.h b/src/image.h index c617330..6b2f4e7 100644 --- a/src/image.h +++ b/src/image.h @@ -42,9 +42,9 @@ typedef struct { typedef struct imagerecord_struct { uint16_t *image; - double tilt; /* Degrees. Defines where the pattern lies in reciprocal space */ - double omega; /* Degrees. Defines where the pattern lies in reciprocal space */ - double slop; /* Degrees. Defines where the pattern lies "on the negative" */ + double tilt; /* Radians. Defines where the pattern lies in reciprocal space */ + double omega; /* Radians. Defines where the pattern lies in reciprocal space */ + double slop; /* Radians. Defines where the pattern lies "on the negative" */ FormulationMode fmode; double pixel_size; diff --git a/src/imagedisplay.c b/src/imagedisplay.c index bc2fa53..00b7285 100644 --- a/src/imagedisplay.c +++ b/src/imagedisplay.c @@ -192,13 +192,13 @@ static gboolean imagedisplay_redraw(GtkWidget *drawingarea, GdkEventExpose *even imagedisplay->imagerecord.y_centre * scale, (imagedisplay->imagerecord.x_centre + imagedisplay->imagerecord.width) * scale, (imagedisplay->imagerecord.y_centre + imagedisplay->imagerecord.width - * tan(deg2rad(imagedisplay->imagerecord.omega))) * scale); + * tan(imagedisplay->imagerecord.omega)) * scale); imagedisplay_draw_line(imagedisplay->gc_tiltaxis, imagedisplay->imagerecord.x_centre * scale, imagedisplay->imagerecord.y_centre * scale, (imagedisplay->imagerecord.x_centre - imagedisplay->imagerecord.width) * scale, (imagedisplay->imagerecord.y_centre - imagedisplay->imagerecord.width - * tan(deg2rad(imagedisplay->imagerecord.omega))) * scale); + * tan(imagedisplay->imagerecord.omega)) * scale); } /* NB This calls the function above, which sorts out stuff */ diff --git a/src/iprtest.c b/src/iprtest.c index 4a06a4a..976ad51 100644 --- a/src/iprtest.c +++ b/src/iprtest.c @@ -111,27 +111,6 @@ int main(int argc, char *argv[]) { GtkWidget *window; GtkWidget *drawing; IPRTestContext *ctx; - double dlx, dly, dlz; - - mapping_rotate(1.0, 0.0, 0.0, &dlx, &dly, &dlz, 0.0, 0.0); - printf("omega=0.0, tilt=0.0 => image 'x' = { %f %f %f }\n", dlx, dly, dlz); - mapping_rotate(0.0, 1.0, 0.0, &dlx, &dly, &dlz, 0.0, 0.0); - printf("omega=0.0, tilt=0.0 => image 'y' = { %f %f %f }\n", dlx, dly, dlz); - - mapping_rotate(1.0, 0.0, 0.0, &dlx, &dly, &dlz, deg2rad(90.0), 0.0); - printf("omega=90.0, tilt=0.0 => image 'x' = { %f %f %f }\n", dlx, dly, dlz); - mapping_rotate(0.0, 1.0, 0.0, &dlx, &dly, &dlz, deg2rad(90.0), 0.0); - printf("omega=90.0, tilt=0.0 => image 'y' = { %f %f %f }\n", dlx, dly, dlz); - - mapping_rotate(1.0, 0.0, 0.0, &dlx, &dly, &dlz, 0.0, deg2rad(90.0)); - printf("omega=0.0, tilt=90.0 => image 'x' = { %f %f %f }\n", dlx, dly, dlz); - mapping_rotate(0.0, 1.0, 0.0, &dlx, &dly, &dlz, 0.0, deg2rad(90.0)); - printf("omega=0.0, tilt=90.0 => image 'y' = { %f %f %f }\n", dlx, dly, dlz); - - mapping_rotate(1.0, 0.0, 0.0, &dlx, &dly, &dlz, deg2rad(90.0), deg2rad(90.0)); - printf("omega=90.0, tilt=90.0 => image 'x' = { %f %f %f }\n", dlx, dly, dlz); - mapping_rotate(0.0, 1.0, 0.0, &dlx, &dly, &dlz, deg2rad(90.0), deg2rad(90.0)); - printf("omega=90.0, tilt=90.0 => image 'y' = { %f %f %f }\n", dlx, dly, dlz); gtk_init(&argc, &argv); @@ -139,8 +118,8 @@ int main(int argc, char *argv[]) { ctx->image = malloc(sizeof(ImageRecord)); ctx->image->image = NULL; - ctx->image->tilt = deg2rad(90.0); - ctx->image->omega = deg2rad(90.0); + ctx->image->tilt = deg2rad(45.0); + ctx->image->omega = deg2rad(0.0); ctx->image->slop = deg2rad(0.0); ctx->image->fmode = FORMULATION_PIXELSIZE; ctx->image->pixel_size = 0.6e8; @@ -156,15 +135,15 @@ int main(int argc, char *argv[]) { ctx->cell = malloc(sizeof(Basis)); /* The "true" cell */ - ctx->cell->a.x = 5.0e9; ctx->cell->a.y = 0.0; ctx->cell->a.z = 0.0; - ctx->cell->b.x = 0.0; ctx->cell->b.y = 5.0e9; ctx->cell->b.z = 0.0; - ctx->cell->c.x = 0.0; ctx->cell->c.y = 0.0; ctx->cell->c.z = 5.0e9; + ctx->cell->a.x = 5.0e9; ctx->cell->a.y = 0.0e9; ctx->cell->a.z = 0.0e9; + ctx->cell->b.x = 0.0e9; ctx->cell->b.y = 5.0e9; ctx->cell->b.z = 0.0e9; + ctx->cell->c.x = 0.0e9; ctx->cell->c.y = 0.0e9; ctx->cell->c.z = 5.0e9; ctx->reflections = reflection_list_from_cell(ctx->cell); ctx->image->features = reproject_get_reflections(ctx->image, ctx->reflections); /* The "model" cell */ - ctx->cell->a.x = 5.2e9; ctx->cell->a.y = 0.1e9; ctx->cell->a.z = 0.0; - ctx->cell->b.x = 0.1e9; ctx->cell->b.y = 4.5e9; ctx->cell->b.z = 0.0; - ctx->cell->c.x = 0.0; ctx->cell->c.y = 0.0; ctx->cell->c.z = 5.0e9; + ctx->cell->a.x = 5.2e9; ctx->cell->a.y = 0.1e9; ctx->cell->a.z = 0.0e9; + ctx->cell->b.x = 0.2e9; ctx->cell->b.y = 4.8e9; ctx->cell->b.z = 0.0e9; + ctx->cell->c.x = 0.0e9; ctx->cell->c.y = 0.0e9; ctx->cell->c.z = 5.0e9; ctx->gen = reflection_list_from_cell(ctx->cell); ctx->image->rflist = reproject_get_reflections(ctx->image, ctx->gen); diff --git a/src/mapping.c b/src/mapping.c index bc6e062..c9872f9 100644 --- a/src/mapping.c +++ b/src/mapping.c @@ -19,6 +19,7 @@ #include "image.h" #include "displaywindow.h" #include "cache.h" +#include "utils.h" void mapping_rotate(double x, double y, double z, double *ddx, double *ddy, double *ddz, double omega, double tilt) { @@ -66,8 +67,8 @@ int mapping_map_to_space(ImageFeature *refl, double *ddx, double *ddy, double *d imagerecord = refl->parent; x = refl->x - imagerecord->x_centre; y = refl->y - imagerecord->y_centre; - tilt = 2*M_PI*(imagerecord->tilt/360); /* Convert to Radians */ - omega = 2*M_PI*(imagerecord->omega/360); /* Likewise */ + tilt = imagerecord->tilt; + omega = imagerecord->omega; k = 1/imagerecord->lambda; /* Calculate an angular description of the reflection */ @@ -162,7 +163,8 @@ void mapping_adjust_axis(ControlContext *ctx, double offset) { int i; for ( i=0; iimages->n_images; i++ ) { - printf("Image #%3i: old omega=%f deg, new omega=%f deg\n", i, ctx->images->images[i].omega, ctx->images->images[i].omega+offset); + printf("Image #%3i: old omega=%f deg, new omega=%f deg\n", i, rad2deg(ctx->images->images[i].omega), + rad2deg(ctx->images->images[i].omega+offset)); ctx->images->images[i].omega += offset; } diff --git a/src/mrc.c b/src/mrc.c index b68166f..57df51c 100644 --- a/src/mrc.c +++ b/src/mrc.c @@ -79,7 +79,7 @@ int mrc_read(ControlContext *ctx) { } else { ctx->lambda = lambda(1000*ext[i].voltage); } - ctx->omega = ext[i].tilt_axis; + ctx->omega = deg2rad(ext[i].tilt_axis); ctx->pixel_size = ext[i].pixel_size; fseek(fh, mrc.next + sizeof(MRCHeader) + mrc.nx*mrc.ny*2*i, SEEK_SET); @@ -93,7 +93,7 @@ int mrc_read(ControlContext *ctx) { } free(image); - image_add(ctx->images, uimage, mrc.nx, mrc.ny, ext[i].a_tilt, ctx); + image_add(ctx->images, uimage, mrc.nx, mrc.ny, deg2rad(ext[i].a_tilt), ctx); } diff --git a/src/qdrp.c b/src/qdrp.c index 5c7df98..aae0c2e 100644 --- a/src/qdrp.c +++ b/src/qdrp.c @@ -20,6 +20,7 @@ #include "readpng.h" #include "control.h" #include "reflections.h" +#include "utils.h" static void qdrp_chomp(char *str) { @@ -65,9 +66,9 @@ static int qdrp_parsefileline(ControlContext *ctx, const char *line) { skip_chars; file = strndup(line+mark, pos-mark); - tilt = strtod(tilt_s, NULL); + tilt = deg2rad(strtod(tilt_s, NULL)); free(tilt_s); - printf("QD: Reading file: Tilt = %f deg, File='%s'\n", tilt, file); + printf("QD: Reading file: Tilt = %f deg, File='%s'\n", rad2deg(tilt), file); if ( readpng_read(file, tilt, ctx) ) { printf("Reconstruction failed.\n"); @@ -102,7 +103,7 @@ static int qdrp_parseline(ControlContext *ctx, const char *line) { } if ( !ctx->omega_set ) { - printf("QD: Parameter 'omega' not specified, and not using tilt-axis-search.\n"); + printf("QD: Parameter 'omega' not specified.\n"); return 1; } @@ -175,10 +176,10 @@ static int qdrp_parseline(ControlContext *ctx, const char *line) { mark = pos; skip_chars; omega_s = strndup(line+mark, pos-mark); - ctx->omega = strtod(omega_s, NULL); + ctx->omega = deg2rad(strtod(omega_s, NULL)); free(omega_s); ctx->omega_set = 1; - printf("QD: omega = %f deg\n", ctx->omega); + printf("QD: omega = %f deg\n", rad2deg(ctx->omega)); } diff --git a/src/refine.c b/src/refine.c index 3a01a9d..fe932d0 100644 --- a/src/refine.c +++ b/src/refine.c @@ -192,6 +192,11 @@ ImageFeature *refine_fit_image(Basis *cell, ImageRecord *image) { /* Prepare for solving */ refine_fix_unconstrained(M); + gsl_matrix_set(M, 2, 0, 0.0); + gsl_matrix_set(M, 2, 1, 0.0); + gsl_matrix_set(M, 2, 2, 1.0); + gsl_matrix_set(M, 0, 2, 0.0); + gsl_matrix_set(M, 1, 2, 0.0); matrix_vector_show(M, p, "RF: "); /* Calculate and display the determinant */ diff --git a/src/reproject.c b/src/reproject.c index 41ec8c5..4036863 100644 --- a/src/reproject.c +++ b/src/reproject.c @@ -58,8 +58,8 @@ ImageFeatureList *reproject_get_reflections(ImageRecord *image, ReflectionList * flist = image_feature_list_new(); - tilt = deg2rad(image->tilt); - omega = deg2rad(image->omega); + tilt = image->tilt; + omega = image->omega; /* Calculate the (normalised) incident electron wavevector * n is rotated "into" the reconstruction, so only one omega step. */ -- cgit v1.2.3