diff options
Diffstat (limited to 'src/luzzatti.c')
-rw-r--r-- | src/luzzatti.c | 91 |
1 files changed, 91 insertions, 0 deletions
diff --git a/src/luzzatti.c b/src/luzzatti.c new file mode 100644 index 0000000..4aca763 --- /dev/null +++ b/src/luzzatti.c @@ -0,0 +1,91 @@ +/* + * luzzatti.c + * + * Plot of R-factor against resolution + * + * (c) 2006-2007 Thomas White <taw27@cam.ac.uk> + * + * synth2d - Two-Dimensional Crystallographic Fourier Synthesis + * + */ + +#ifdef HAVE_CONFIG_H +#include <config.h> +#endif + +#include <stdlib.h> +#include <stdio.h> +#include <math.h> +#include <string.h> + +#include "data.h" +#include "main.h" +#include "statistics.h" +#include "displaywindow.h" +#include "model.h" + +#define LUZZATTI_BINS 50 + +static void luzzatti_gpwrite(FILE *gnuplot, const char *string) { + fwrite(string, strlen(string), 1, gnuplot); +} + +void luzzatti_show() { + + FILE *fh; + float a = data_a(); + float b = data_b(); + float c = data_c(); + unsigned int i; + FILE *gnuplot; + double scale; + ReflectionList *reflections; + ReflectionList *model_reflections; + + reflections = main_reflist(); + model_reflections = model_calculate_f(reflections, NULL, 69); + + scale = stat_scale(reflections, model_reflections); + fh = fopen("synth2d-luzzatti.dat", "w"); + + for ( i=1; i<reflections->n_reflections; i++ ) { + + double s, residual; + signed int h, k, l; + + h = reflections->refs[i].h; + k = reflections->refs[i].k; + l = reflections->refs[i].l; + + s = sqrt(((h*h)/(a*a)) + ((k*k)/(b*b)) + ((l*l)/(c*c))); + residual = 100*fabs((reflections->refs[i].amplitude - scale*model_reflections->refs[i].amplitude)/reflections->refs[i].amplitude); + printf("%3i %3i %3i: obs=%f calc=%f => residual=%f\n", h, k, l, reflections->refs[i].amplitude, + scale*model_reflections->refs[i].amplitude, residual); + fprintf(fh, "%f %f\n", s, residual); + } + + fclose(fh); + + gnuplot = popen("gnuplot -persist -", "w"); + if ( !gnuplot ) { + error_report("Couldn't invoke gnuplot. Please check your PATH."); + return; + } + + luzzatti_gpwrite(gnuplot, "set autoscale\n"); + luzzatti_gpwrite(gnuplot, "unset log\n"); + luzzatti_gpwrite(gnuplot, "unset label\n"); + luzzatti_gpwrite(gnuplot, "set xtic auto\n"); + luzzatti_gpwrite(gnuplot, "set ytic auto\n"); + luzzatti_gpwrite(gnuplot, "set grid\n"); + luzzatti_gpwrite(gnuplot, "set ylabel 'R (%)' font \"Helvetica,10\"\n"); + luzzatti_gpwrite(gnuplot, "set xlabel '1/d / nm^(-1)' font \"Helvetica,10\"\n"); + luzzatti_gpwrite(gnuplot, "set title 'Luzzatti Plot' font \"Helvetica,10\"\n"); + luzzatti_gpwrite(gnuplot, "plot 'synth2d-luzzatti.dat'\n"); + + if ( pclose(gnuplot) == -1 ) { + error_report("gnuplot returned an error code."); + return; + } + +} |