/*
 * amplitude-r.c
 *
 * Plot of R-factor against amplitude
 *
 * (c) 2006 Thomas White <taw27@cam.ac.uk>
 *  Synth2d - two-dimensional 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"

static void amplituder_gpwrite(FILE *gnuplot, const char *string) {
	fwrite(string, strlen(string), 1, gnuplot);
}

void amplituder_show() {

	FILE *fh;
	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-amplitude-r.dat", "w");
		
	for ( i=1; i<reflections->n_reflections; i++ ) {

		double residual;
		signed int h, k, l;
		
		h = reflections->refs[i].h;
		k = reflections->refs[i].k;
		l = reflections->refs[i].l;
		
		residual = 100*fabs((reflections->refs[i].amplitude - scale*model_reflections->refs[i].amplitude)/reflections->refs[i].amplitude);
		if ( !isinf(residual) ) {
			fprintf(fh, "%f %f\n", reflections->refs[i].amplitude, residual);
		}

	}
	
	fclose(fh);
	
	gnuplot = popen("gnuplot -persist -", "w");
	if ( !gnuplot ) {
		error_report("Couldn't invoke gnuplot.  Please check your PATH.");
		return;
	}
	
	amplituder_gpwrite(gnuplot, "set autoscale\n");
	amplituder_gpwrite(gnuplot, "unset log\n");
	amplituder_gpwrite(gnuplot, "unset label\n");
	amplituder_gpwrite(gnuplot, "set xtic auto\n");
	amplituder_gpwrite(gnuplot, "set ytic auto\n");
	amplituder_gpwrite(gnuplot, "set grid\n");
	amplituder_gpwrite(gnuplot, "set ylabel 'R (%)' font \"Helvetica,10\"\n");
	amplituder_gpwrite(gnuplot, "set xlabel '|Fobs|' font \"Helvetica,10\"\n");
	amplituder_gpwrite(gnuplot, "set title 'Amplitude-R Plot' font \"Helvetica,10\"\n");
	amplituder_gpwrite(gnuplot, "plot 'synth2d-amplitude-r.dat'\n");
	amplituder_gpwrite(gnuplot, "replot\n");
	
	if ( pclose(gnuplot) == -1 ) {
		error_report("gnuplot returned an error code.");
		return;
	}
	
}