/*	**************************************************************	*/

/*
 *	Module CALIB
 *	This module contains utility functions to apply calibration
 *	data to measurements.
 *
 *	Zachary Wolf
 *	12/13/02
 */

/*	**************************************************************	*/

/*	INCLUDES  */
#include <analysis.h>
#include <userint.h>
#include <formatio.h>
#include <utility.h>
#include <ansi_c.h>
#include "calib.h"

/*	**************************************************************	*/

/*	PRIVATE GLOBAL VARIABLES  */
static char msg[CALIB_MAX_NUM_CHAR];

/*	**************************************************************	*/

/*	PRIVATE FUNCTIONS  */
void calib_error(char* message);
void calib_message(char* message);
void calib_sort_calib_data(struct calib_data_struct* calib_data);

/*	**************************************************************	*/

/*	PUBLIC FUNCTIONS  */

/*	**************************************************************	*/

/*	CALIB DATA  */

/*	**************************************************************	*/
/*
 *	calib_get_calib_data
 *	This function opens a file and reads the calibration data
 *	from the file.  The data is returned in the calibration data
 *	structure.
 *
 *	Input:
 *	calib_data_file_name, name of the calibration file
 *
 *	Output:
 *	calib_data, structure containing the calibration data
 *
 *	Zachary Wolf
 *	12/13/02
 */

void calib_get_calib_data(char calib_data_file_name[], struct calib_data_struct* calib_data)
{
	/* Declare variables */
	FILE* file_ptr;
	char buf[CALIB_MAX_NUM_CHAR];
	int num_calib_pts;
	double meas_calib, real_calib;
	int err;
	int i;
	
	/* Open the calibration file */
	file_ptr = fopen (calib_data_file_name, "r");
	if (file_ptr == NULL)
	{
		calib_error("Could not open calibration file.");
		return;
	}

	/* Name of the probe */
	fgets(buf, CALIB_MAX_NUM_CHAR, file_ptr);
	if (strstr(buf, "probe_name") == NULL)
	{
		calib_error("Error reading name of the probe!");
		return;
	}
	err = fscanf(file_ptr, "%s\n", buf);
	if (err != 1)
	{
		calib_error("Error reading probe name!");
		return;
	}
	strcpy(calib_data->probe_name, buf);
	
	/* Number of calibration points */
	fgets(buf, CALIB_MAX_NUM_CHAR, file_ptr);
	if (strstr(buf, "num_calib_pts") == NULL)
	{
		calib_error("Error reading num_calib_points!");
		return;
	}
	err = fscanf(file_ptr, "%i\n", &num_calib_pts);
	if (err != 1)
	{
		calib_error("Error reading num_calib_pts!");
		return;
	}
	calib_data->num_calib_pts = num_calib_pts;
	
	/* Check */
	if (num_calib_pts > CALIB_MAX_NUM_CALIB_PTS)
	{
		calib_error("Error, num_calib_pts is too large");
		return;
	}
	
	/* Calibration data */
	fgets(buf, CALIB_MAX_NUM_CHAR, file_ptr);
	if (strstr(buf, "meas_real_calib_data") == NULL)
	{
		calib_error("Error reading calib_data!");
		return;
	}
	for (i = 0; i < num_calib_pts; i++)
	{
		err = fscanf(file_ptr, "%lf %lf\n", &meas_calib, &real_calib);
		if (err != 2)
		{
			calib_error("Error reading calibration data!");
			return;
		}
		calib_data->meas_calib[i] = meas_calib;
		calib_data->real_calib[i] = real_calib;
	}
	
	/* Close the calibration file */
	fclose(file_ptr);
	
	/* Sort the calibration data from smallest meas_calib to largest */
	calib_sort_calib_data(calib_data);
	
	/* Done */
	return;
}

/*	**************************************************************	*/
/*
 *	calib_plot_calib_data
 *	This function plots the calibration data points.
 *
 *	Input:
 *	calib_data, structure containing the calibration data
 *
 *	Zachary Wolf
 *	12/13/02
 */

void calib_plot_calib_data(struct calib_data_struct calib_data)
{
	/* Declare variables */
	FILE* file_ptr;
	int i;
	
	/* Open the plt file */
	file_ptr = fopen("Val_calib_data.plt", "w");
	if (file_ptr == NULL) return;
	
	/* Write the values to the plt file */
	fprintf(file_ptr, "%% meas_calib, real_calib\n");
	for (i = 0; i < calib_data.num_calib_pts; i++)
		fprintf(file_ptr, "%f  %f\n", calib_data.meas_calib[i], calib_data.real_calib[i]);
	
	/* Close the plt file */
	fclose(file_ptr);
	
	/* Done */
	return;
}

/*	**************************************************************	*/
/*
 *	calib_write_calib_data
 *	This function writes the calibration data to a file.
 *
 *	Input:
 *	calib_data_file_name, name of the calibration file to write
 *	calib_data, structure containing the calibration data
 *
 *	Zachary Wolf
 *	12/13/02
 */

void calib_write_calib_data(char calib_data_file_name[], struct calib_data_struct calib_data)
{
	/* Declare variables */
	FILE* file_ptr;
	int i;
	
	/* Create the calibration file */
	file_ptr = fopen (calib_data_file_name, "w");
	if (file_ptr == NULL)
	{
		calib_error("Could not create calibration file.");
		return;
	}

	/* Name of the probe */
	fprintf(file_ptr, "probe_name\n");
	fprintf(file_ptr, "%s\n", calib_data.probe_name);
	
	/* Number of calibration points */
	fprintf(file_ptr, "num_calib_pts\n");
	fprintf(file_ptr, "%i\n", calib_data.num_calib_pts);
	
	/* Calibration fit coefficients */
	fprintf(file_ptr, "meas_real_calib_data\n");
	for (i = 0; i < calib_data.num_calib_pts; i++)
	{
		fprintf(file_ptr, "%f %f\n", calib_data.meas_calib[i], calib_data.real_calib[i]);
	}
	
	/* Close the calibration file */
	fclose(file_ptr);
	
	/* Done */
	return;
}

/*	**************************************************************	*/

/*	CUBIC SPLINE INTERPOLATION  */

/*	**************************************************************	*/
/*
 *	calib_apply_calib_data_cubic_spline_interp
 *	This function calculates the output value from the measured
 *	value and the calibration data.
 *
 *	Input:
 *	calib_data, structure containing the calibration data
 *	num_meas, number of measurements
 *	meas_val[0 to num_meas - 1], measured values (eg. V)
 *
 *	Output:
 *	calib_out[0 to num_meas - 1], calibrated output value (eg. T)
 *
 *	Zachary Wolf
 *	12/13/02
 */

void calib_apply_calib_data_cubic_spline_interp(struct calib_data_struct calib_data, int num_meas, double meas_val[], double calib_out[])
{
	/* Declare variables */
	int i;
	double x[4], y[4], z[4], coef[4], mse;
	int n;
	double b1, b2;
	double y2[CALIB_MAX_NUM_CALIB_PTS];
	int status;

	/* Check that the calibration data has been read */
	if (strcmp(calib_data.probe_name, "") == 0) calib_error("No calibration data file has been read");
	
	/* Make sure each measured value is within the range of the calibration */
	for (i = 0; i < num_meas; i++)
	{
		if (meas_val[i] < calib_data.meas_calib[0] || meas_val[i] > calib_data.meas_calib[calib_data.num_calib_pts - 1])
			calib_error("Insufficient calibration data");
	}
	
	/* Handle the ends of the calibration range */
	/* Estimate the first derivative at the small meas_calib end */
	for (i = 0; i < 2; i++)
	{
		x[i] = calib_data.meas_calib[i];
		y[i] = calib_data.real_calib[i];   
	}
	status = PolyFit(x, y, 2, 1, z, coef, &mse);
	if (status != 0) calib_error("Problem with the cubic spline");
	b1 = coef[1];
	/* Estimate the first derivative at the large meas_calib end */
	for (i = 0; i < 2; i++)
	{
		x[i] = calib_data.meas_calib[calib_data.num_calib_pts - 1 - 1 + i];
		y[i] = calib_data.real_calib[calib_data.num_calib_pts - 1 - 1 + i];
	}
	status = PolyFit(x, y, 2, 1, z, coef, &mse);
	if (status != 0) calib_error("Problem with the cubic spline");
	b2 = coef[1];
	
	/* Get coefficients for a cubic spline interpolation */
	status = Spline(calib_data.meas_calib, calib_data.real_calib, calib_data.num_calib_pts, b1, b2, y2);
	if (status != 0) calib_error("Problem with the cubic spline");
	
	/* Do the interpolation for each measurement */
	for (i = 0; i < num_meas; i++)
	{
		status = SpInterp(calib_data.meas_calib, calib_data.real_calib, y2, calib_data.num_calib_pts, meas_val[i], &calib_out[i]);
		if (status != 0) calib_error("Problem with the cubic spline");
	}
	
	/* Done */
	return;
}

/*	**************************************************************	*/
/*
 *	calib_plot_cubic_spline_interp
 *	This function plots the results of the cubic spline interpolation
 *	between the calibration points.
 *
 *	Input:
 *	calib_data, structure containing the calibration data
 *
 *	Zachary Wolf
 *	12/13/02
 */

void calib_plot_cubic_spline_interp(struct calib_data_struct calib_data)
{
	/* Declare variables */
	int num_meas = 201;
	double meas_val[201];
	double calib_out[201];
	double range_min, range_max;
	double delta_meas_val;
	int i;
	FILE* file_ptr;
	
	/* Generate a number of points within the range of the calibration data */
	range_min = calib_data.meas_calib[0];
	range_max = calib_data.meas_calib[calib_data.num_calib_pts - 1];
	delta_meas_val = (range_max - range_min) / (num_meas - 1);
	for (i = 0; i < num_meas; i++) meas_val[i] = range_min + i * delta_meas_val; 
	
	/* Apply the calibration function to the generated measured points */
	calib_apply_calib_data_cubic_spline_interp(calib_data, num_meas, meas_val, calib_out);
	
	/* Plot the results */
	/* Open the plt file */
	file_ptr = fopen("Val_cubic_spline_interp.plt", "w");
	if (file_ptr == NULL) return;
	
	/* Write the values to the plt file */
	fprintf(file_ptr, "%% meas_val, calib_out\n");
	for (i = 0; i < num_meas; i++)
		fprintf(file_ptr, "%f  %f\n", meas_val[i], calib_out[i]);
	
	/* Close the plt file */
	fclose(file_ptr);
	
	/* Done */
	return;
}

/*	**************************************************************	*/
/*
 *	calib_test_cubic_spline_interp
 *	This function performs tests of the cubic spline interpolation functions.
 *
 *	Zachary Wolf
 *	12/13/02
 */

void calib_test_cubic_spline_interp(void)
{
	/* Declare variables */
	struct calib_data_struct gen_calib_data;
	struct calib_data_struct calib_data;
	struct calib_data_struct meas_data;
	int i;
	double calib_out[CALIB_MAX_NUM_CALIB_PTS];
	double error[CALIB_MAX_NUM_CALIB_PTS];
	FILE* file_ptr;
	
	/* Generate fake calibration data */
	calib_generate_test_calib_data(20, -2., 2., &gen_calib_data);
	
	/* Write the fake calibration data to a file */
	calib_write_calib_data("Val_test_calib_data_file.dat", gen_calib_data);
	
	/* Read the fake calibration data from the file */
	calib_get_calib_data("Val_test_calib_data_file.dat", &calib_data);
	
	/* Check the probe name */
	if (strcmp(calib_data.probe_name, "test_calib") != 0) printf("Problem with the calibration data.\n");
	
	/* Plot the fake calibration data */
	calib_plot_calib_data(calib_data);
	
	/* Plot the cubic spline interpolation of the calibration data */
	calib_plot_cubic_spline_interp(calib_data);
	
	/* Generate fake measurement data */
	calib_generate_test_calib_data(100, -2., 2., &meas_data);
	
	/* Apply the calibration to the fake measurement data */
	calib_apply_calib_data_cubic_spline_interp(calib_data, meas_data.num_calib_pts, meas_data.meas_calib, calib_out);
	
	/* Compare the calibrated output to the expected value */
	for (i = 0; i < meas_data.num_calib_pts; i++)
	{
		error[i] = meas_data.real_calib[i] - calib_out[i];
	}
	
	/* Record the results */
	file_ptr = fopen("Val_test_cubic_spline_interp.plt", "w");
	if (file_ptr == NULL) return;
	
	/* Write the values to the plt file */
	fprintf(file_ptr, "%% meas_val, true_val, calib_out, error\n");
	for (i = 0; i < meas_data.num_calib_pts; i++)
		fprintf(file_ptr, "%f  %f  %f  %f\n", meas_data.meas_calib[i], meas_data.real_calib[i], calib_out[i], error[i]);
	
	/* Close the plt file */
	fclose(file_ptr);
	
	/* Done */
	return;
}

/*	**************************************************************	*/

/*	POLYNOMIAL FIT INTERPOLATION  */

/*	**************************************************************	*/
/*
 *	calib_get_calib_poly_fit
 *	This function opens a file and reads the parameters of the polynomial fit
 *	to the calibration data from the file.  The data is returned in the calibration
 *	polynomial fit structure.
 *
 *	Input:
 *	calib_poly_fit_file_name, name of the calibration polynomial fit parameter file
 *
 *	Output:
 *	calib_poly_fit, structure containing the calibration polynomial fit parameters
 *
 *	Zachary Wolf
 *	12/13/02
 */

void calib_get_calib_poly_fit(char calib_poly_fit_file_name[], struct calib_poly_fit_struct* calib_poly_fit)
{
	/* Declare variables */
	FILE* file_ptr;
	char buf[CALIB_MAX_NUM_CHAR];
	double meas_range_min, meas_range_max;
	int fit_order;
	double fit_coef;
	int err;
	int i;
	
	/* Open the calibration file */
	file_ptr = fopen (calib_poly_fit_file_name, "r");
	if (file_ptr == NULL)
	{
		calib_error("Could not open calibration file.");
		return;
	}

	/* Name of the probe */
	fgets(buf, CALIB_MAX_NUM_CHAR, file_ptr);
	if (strstr(buf, "probe_name") == NULL)
	{
		calib_error("Error reading name of the probe!");
		return;
	}
	err = fscanf(file_ptr, "%s\n", buf);
	if (err != 1)
	{
		calib_error("Error reading probe name!");
		return;
	}
	strcpy(calib_poly_fit->probe_name, buf);
	
	/* Range min */
	fgets(buf, CALIB_MAX_NUM_CHAR, file_ptr);
	if (strstr(buf, "meas_range_min") == NULL)
	{
		calib_error("Error reading meas_range_min!");
		return;
	}
	err = fscanf(file_ptr, "%lf\n", &meas_range_min);
	if (err != 1)
	{
		calib_error("Error reading meas_range_min!");
		return;
	}
	calib_poly_fit->meas_range_min = meas_range_min;
	
	/* Range max */
	fgets(buf, CALIB_MAX_NUM_CHAR, file_ptr);
	if (strstr(buf, "meas_range_max") == NULL)
	{
		calib_error("Error reading meas_range_max!");
		return;
	}
	err = fscanf(file_ptr, "%lf\n", &meas_range_max);
	if (err != 1)
	{
		calib_error("Error reading meas_range_max!");
		return;
	}
	calib_poly_fit->meas_range_max = meas_range_max;
	
	/* Order of the polynomial fit */
	fgets(buf, CALIB_MAX_NUM_CHAR, file_ptr);
	if (strstr(buf, "fit_order") == NULL)
	{
		calib_error("Error reading fit_order!");
		return;
	}
	err = fscanf(file_ptr, "%i\n", &fit_order);
	if (err != 1)
	{
		calib_error("Error reading fit_order!");
		return;
	}
	calib_poly_fit->fit_order = fit_order;
	
	/* Check */
	if (fit_order > CALIB_MAX_FIT_ORDER_PLUS_1 - 1)
	{
		calib_error("Error, fit_order is too large");
		return;
	}
	
	/* Calibration fit coefficients */
	fgets(buf, CALIB_MAX_NUM_CHAR, file_ptr);
	if (strstr(buf, "fit_coef") == NULL)
	{
		calib_error("Error reading calib_poly_fit_coef!");
		return;
	}
	for (i = 0; i <= fit_order; i++)
	{
		err = fscanf(file_ptr, "%lf\n", &fit_coef);
		if (err != 1)
		{
			calib_error("Error reading calibration fit data!");
			return;
		}
		calib_poly_fit->fit_coef[i] = fit_coef;
	}
	
	/* Close the calibration file */
	fclose(file_ptr);
	
	/* Done */
	return;
}

/*	**************************************************************	*/
/*
 *	calib_apply_calib_poly_fit_interp
 *	This function calculates the output value from the measured
 *	value and the calibration polynomial fit parameters.
 *
 *	Input:
 *	calib_poly_fit, structure containing the calibration polynomial fit parameters
 *	num_meas, number of measurements
 *	meas_val[0 to num_meas - 1], measured values (eg. V)
 *
 *	Output:
 *	calib_out[0 to num_meas - 1], calibrated output value (eg. T)
 *
 *	Zachary Wolf
 *	12/13/02
 */

void calib_apply_calib_poly_fit_interp(struct calib_poly_fit_struct calib_poly_fit, int num_meas, double meas_val[], double calib_out[])
{
	/* Declare variables */
	int i;
	int status;

	/* Check that the calibration fit data has been read */
	if (strcmp(calib_poly_fit.probe_name, "") == 0) calib_error("No calibration fit file has been read");
	
	/* Make sure each measured value is within the range of the calibration */
	for (i = 0; i < num_meas; i++)
	{
		if (meas_val[i] < calib_poly_fit.meas_range_min || meas_val[i] > calib_poly_fit.meas_range_max)
			calib_error("Insufficient calibration data");
	}
	
	/* Apply the fit to the measured values */
	status = PolyEv1D(meas_val, num_meas, calib_poly_fit.fit_coef, calib_poly_fit.fit_order + 1, calib_out);
	if (status != 0) calib_error("Problem with the polynomial fit");
	
	/* Done */
	return;
}

/*	**************************************************************	*/
/*
 *	calib_plot_poly_fit_interp
 *	This function plots the results of the polynomial fit interpolation
 *	between the calibration points.
 *
 *	Input:
 *	calib_poly_fit, structure containing the polynomial fit coefficients to the calibration data
 *
 *	Zachary Wolf
 *	12/13/02
 */

void calib_plot_poly_fit_interp(struct calib_poly_fit_struct calib_poly_fit)
{
	/* Declare variables */
	int num_meas = 201;
	double meas_val[201];
	double calib_out[201];
	double range_min, range_max;
	double delta_meas_val;
	int i;
	FILE* file_ptr;
	
	/* Generate a number of points within the range of the calibration data */
	range_min = calib_poly_fit.meas_range_min;
	range_max = calib_poly_fit.meas_range_max;
	delta_meas_val = (range_max - range_min) / (num_meas - 1);
	for (i = 0; i < num_meas; i++) meas_val[i] = range_min + i * delta_meas_val; 
	
	/* Apply the calibration function to the generated measured points */
	calib_apply_calib_poly_fit_interp(calib_poly_fit, num_meas, meas_val, calib_out);
	
	/* Plot the results */
	/* Open the plt file */
	file_ptr = fopen("Val_poly_fit_interp.plt", "w");
	if (file_ptr == NULL) return;
	
	/* Write the values to the plt file */
	fprintf(file_ptr, "%% meas_val, calib_out\n");
	for (i = 0; i < num_meas; i++)
		fprintf(file_ptr, "%f  %f\n", meas_val[i], calib_out[i]);
	
	/* Close the plt file */
	fclose(file_ptr);
	
	/* Done */
	return;
}

/*	**************************************************************	*/
/*
 *	calib_calc_poly_poly_fit_coef
 *	This function calculates the  coefficients in a polynomial fit
 *	to the calibration data.
 *
 *	Input:
 *	calib_data, structure containing the calibration data
 *
 *	Output:
 *	calib_poly_fit, structure containing the polynomial fit coefficients
 *
 *	Zachary Wolf
 *	12/13/02
 */

void calib_calc_poly_fit_coef(struct calib_data_struct calib_data, struct calib_poly_fit_struct* calib_poly_fit)
{
	/* Declare variables */
	int fit_order = 6;
	double z[CALIB_MAX_NUM_CALIB_PTS];
	double mse;
	int status;

	/* Put the probe name in the calib_poly_fit structure */
	strcpy(calib_poly_fit->probe_name, calib_data.probe_name);
	
	/* Fill in the range of applicability of the fit */
	calib_poly_fit->meas_range_min = calib_data.meas_calib[0];
	calib_poly_fit->meas_range_max = calib_data.meas_calib[calib_data.num_calib_pts - 1];
	
	/* Fill in the order of the fit */
	calib_poly_fit->fit_order = fit_order;
	
	/* Fit the calibration data */
	status = PolyFit(calib_data.meas_calib, calib_data.real_calib, calib_data.num_calib_pts, fit_order, z, calib_poly_fit->fit_coef, &mse);
	if (status != 0) calib_error("Problem with the polynomial fit");
	
	/* Done */
	return;
}

/*	**************************************************************	*/
/*
 *	calib_write_poly_fit_coef
 *	This function writes the polynomial fit coefficients to a file.
 *
 *	Input:
 *	calib_poly_fit_file_name, name of the calibration file to write
 *	calib_poly_fit, structure containing the polynomial fit information
 *
 *	Zachary Wolf
 *	12/13/02
 */

void calib_write_poly_fit_coef(char calib_poly_fit_file_name[], struct calib_poly_fit_struct calib_poly_fit)
{
	/* Declare variables */
	FILE* file_ptr;
	int i;
	
	/* Create the calibration file */
	file_ptr = fopen (calib_poly_fit_file_name, "w");
	if (file_ptr == NULL)
	{
		calib_error("Could not create calibration file.");
		return;
	}

	/* Name of the probe */
	fprintf(file_ptr, "probe_name\n");
	fprintf(file_ptr, "%s\n", calib_poly_fit.probe_name);
	
	/* Range min */
	fprintf(file_ptr, "meas_range_min\n");
	fprintf(file_ptr, "%f\n", calib_poly_fit.meas_range_min);
	
	/* Range max */
	fprintf(file_ptr, "meas_range_max\n");
	fprintf(file_ptr, "%f\n", calib_poly_fit.meas_range_max);
	
	/* Order of the polynomial fit */
	fprintf(file_ptr, "fit_order\n");
	fprintf(file_ptr, "%i\n", calib_poly_fit.fit_order);
	
	/* Calibration fit coefficients */
	fprintf(file_ptr, "fit_coef\n");
	for (i = 0; i <= calib_poly_fit.fit_order; i++)
	{
		fprintf(file_ptr, "%f\n", calib_poly_fit.fit_coef[i]);
	}
	
	/* Close the calibration file */
	fclose(file_ptr);
	
	/* Done */
	return;
}

/*	**************************************************************	*/
/*
 *	calib_plot_poly_fit_residuals
 *	This function plots the residuals of the polynomial fit to the calibration
 *	data at the calibration points.
 *
 *	Input:
 *	calib_data, structure containing the calibration data
 *	calib_poly_fit, structure containing the polynomial fit coefficients
 *
 *	Zachary Wolf
 *	12/13/02
 */

void calib_plot_poly_fit_residuals(struct calib_data_struct calib_data, struct calib_poly_fit_struct calib_poly_fit)
{
	/* Declare variables */
	int i;
	double calib_out[CALIB_MAX_NUM_CALIB_PTS];
	double residual[CALIB_MAX_NUM_CALIB_PTS];
	FILE* file_ptr;

	/* Check that the calibration data has been read */
	if (strcmp(calib_data.probe_name, "") == 0) calib_error("No calibration data file has been read");
	if (strcmp(calib_poly_fit.probe_name, "") == 0) calib_error("No calibration fit file has been read");
	
	/* Apply the calibration function to the measured calibration points */
	calib_apply_calib_poly_fit_interp(calib_poly_fit, calib_data.num_calib_pts, calib_data.meas_calib, calib_out);
	
	/* Compare the calibrated output to the expected value */
	for (i = 0; i < calib_data.num_calib_pts; i++)
	{
		residual[i] = calib_data.real_calib[i] - calib_out[i];
	}
	
	/* Plot the results */
	/* Open the plt file */
	file_ptr = fopen("Val_poly_fit_residuals.plt", "w");
	if (file_ptr == NULL) return;
	
	/* Write the values to the plt file */
	fprintf(file_ptr, "%% meas_calib, real_calib, calib_out, residual\n");
	for (i = 0; i < calib_data.num_calib_pts; i++)
		fprintf(file_ptr, "%f  %f  %f  %f\n", calib_data.meas_calib[i], calib_data.real_calib[i], calib_out[i], residual[i]);
	
	/* Close the plt file */
	fclose(file_ptr);
	
	/* Done */
	return;
}

/*	**************************************************************	*/
/*
 *	calib_test_poly_fit_interp
 *	This function performs tests of the polynomial fit interpolation functions.
 *
 *	Zachary Wolf
 *	12/13/02
 */

void calib_test_poly_fit_interp(void)
{
	/* Declare variables */
	struct calib_data_struct gen_calib_data;
	struct calib_data_struct calib_data;
	struct calib_poly_fit_struct gen_calib_poly_fit;
	struct calib_poly_fit_struct calib_poly_fit;
	struct calib_data_struct meas_data;
	int i;
	double calib_out[CALIB_MAX_NUM_CALIB_PTS];
	double error[CALIB_MAX_NUM_CALIB_PTS];
	FILE* file_ptr;
	
	/* Generate fake calibration data */
	calib_generate_test_calib_data(20, -2., 2., &gen_calib_data);
	
	/* Write the fake calibration data to a file */
	calib_write_calib_data("Val_test_calib_data_file.dat", gen_calib_data);
	
	/* Read the fake calibration data from the file */
	calib_get_calib_data("Val_test_calib_data_file.dat", &calib_data);
	
	/* Check the probe name */
	if (strcmp(calib_data.probe_name, "test_calib") != 0) printf("Problem with the calibration data.\n");
	
	/* Plot the fake calibration data */
	calib_plot_calib_data(calib_data);
	
	/* Calculate the polynomial coefficients in a fit to the calibration data */
	calib_calc_poly_fit_coef(calib_data, &gen_calib_poly_fit);
	
	/* Plot the residuals of the polynomial fit */
	calib_plot_poly_fit_residuals(calib_data, gen_calib_poly_fit);

	/* Write the polynomial coefficients to a file */
	calib_write_poly_fit_coef("Val_test_calib_poly_fit_file.dat", gen_calib_poly_fit);
	
	/* Read the polynomial fit coefficients from the file */
	calib_get_calib_poly_fit("Val_test_calib_poly_fit_file.dat", &calib_poly_fit);
	
	/* Plot the polynomial fit interpolation of the calibration data */
	calib_plot_poly_fit_interp(calib_poly_fit);
	
	/* Generate fake measurement data */
	calib_generate_test_calib_data(100, -2., 2., &meas_data);
	
	/* Apply the calibration to the fake measurement data */
	calib_apply_calib_poly_fit_interp(calib_poly_fit, meas_data.num_calib_pts, meas_data.meas_calib, calib_out);
	
	/* Compare the calibrated output to the expected value */
	for (i = 0; i < meas_data.num_calib_pts; i++)
	{
		error[i] = meas_data.real_calib[i] - calib_out[i];
	}
	
	/* Record the results */
	file_ptr = fopen("Val_test_poly_fit_interp.plt", "w");
	if (file_ptr == NULL) return;
	
	/* Write the values to the plt file */
	fprintf(file_ptr, "%% meas_val, true_val, calib_out, error\n");
	for (i = 0; i < meas_data.num_calib_pts; i++)
		fprintf(file_ptr, "%f  %f  %f  %f\n", meas_data.meas_calib[i], meas_data.real_calib[i], calib_out[i], error[i]);
	
	/* Close the plt file */
	fclose(file_ptr);
	
	/* Done */
	return;
}

/*	**************************************************************	*/

/*	STEINHART-HART CALIB */

/*	**************************************************************	*/
/*
 *	calib_get_calib_steinhart_hart
 *	This function opens a file and reads the parameters of the Steinhart-Hart
 *	equation from the file.  The data is returned in the calibration
 *	steinhart_hart structure.  This is used to determine temperature
 *	from resistance measurements.
 *
 *	Input:
 *	calib_steinhart_hart_file_name, name of the calibration parameter file
 *
 *	Output:
 *	calib_steinhart_hart, structure containing the calibration steinhart_hart parameters
 *
 *	Zachary Wolf
 *	12/13/02
 */

void calib_get_calib_steinhart_hart(char calib_steinhart_hart_file_name[], struct calib_steinhart_hart_struct* calib_steinhart_hart)
{
	/* Declare variables */
	FILE* file_ptr;
	char buf[CALIB_MAX_NUM_CHAR];
	double meas_range_min, meas_range_max;
	double a_coef, b_coef, c_coef;
	int err;
	
	/* Open the calibration file */
	file_ptr = fopen (calib_steinhart_hart_file_name, "r");
	if (file_ptr == NULL)
	{
		calib_error("Could not open calibration file.");
		return;
	}

	/* Name of the probe */
	fgets(buf, CALIB_MAX_NUM_CHAR, file_ptr);
	if (strstr(buf, "probe_name") == NULL)
	{
		calib_error("Error reading name of the probe!");
		return;
	}
	err = fscanf(file_ptr, "%s\n", buf);
	if (err != 1)
	{
		calib_error("Error reading probe name!");
		return;
	}
	strcpy(calib_steinhart_hart->probe_name, buf);
	
	/* Range min */
	fgets(buf, CALIB_MAX_NUM_CHAR, file_ptr);
	if (strstr(buf, "meas_range_min") == NULL)
	{
		calib_error("Error reading meas_range_min!");
		return;
	}
	err = fscanf(file_ptr, "%lf\n", &meas_range_min);
	if (err != 1)
	{
		calib_error("Error reading meas_range_min!");
		return;
	}
	calib_steinhart_hart->meas_range_min = meas_range_min;
	
	/* Range max */
	fgets(buf, CALIB_MAX_NUM_CHAR, file_ptr);
	if (strstr(buf, "meas_range_max") == NULL)
	{
		calib_error("Error reading meas_range_max!");
		return;
	}
	err = fscanf(file_ptr, "%lf\n", &meas_range_max);
	if (err != 1)
	{
		calib_error("Error reading meas_range_max!");
		return;
	}
	calib_steinhart_hart->meas_range_max = meas_range_max;
	
	/* Steinhart-Hart coefficients */
	fgets(buf, CALIB_MAX_NUM_CHAR, file_ptr);
	if (strstr(buf, "abc_coef") == NULL)
	{
		calib_error("Error reading abc_coef!");
		return;
	}

	/* A coefficient */
	err = fscanf(file_ptr, "%lf\n", &a_coef);
	if (err != 1)
	{
		calib_error("Error reading calibration data!");
		return;
	}
	calib_steinhart_hart->a_coef = a_coef;

	/* B coefficient */
	err = fscanf(file_ptr, "%lf\n", &b_coef);
	if (err != 1)
	{
		calib_error("Error reading calibration data!");
		return;
	}
	calib_steinhart_hart->b_coef = b_coef;

	/* C coefficient */
	err = fscanf(file_ptr, "%lf\n", &c_coef);
	if (err != 1)
	{
		calib_error("Error reading calibration data!");
		return;
	}
	calib_steinhart_hart->c_coef = c_coef;
	
	/* Close the calibration file */
	fclose(file_ptr);
	
	/* Done */
	return;
}

/*	**************************************************************	*/
/*
 *	calib_apply_calib_steinhart_hart
 *	This function calculates the output temperature from the measured
 *	resistance and the calibration Steinhart-Hart parameters.
 *
 *	Input:
 *	calib_steinhart_hart, structure containing the calibration steinhart_hart parameters
 *	meas_res, measured resistance (ohm)
 *
 *	Output:
 *	temp, calibrated output temperature (deg C)
 *
 *	Zachary Wolf
 *	12/13/02
 */

void calib_apply_calib_steinhart_hart(struct calib_steinhart_hart_struct calib_steinhart_hart, double meas_res, double* temp)
{
	/* Declare variables */
	double a_coef, b_coef, c_coef;
	double R, Tinv, Tk;

	/* Check that the calibration data has been read */
	if (strcmp(calib_steinhart_hart.probe_name, "") == 0) calib_error("No calibration file has been read");
	
	/* Make sure the measured value is within the range of the calibration */
	if (meas_res < calib_steinhart_hart.meas_range_min || meas_res > calib_steinhart_hart.meas_range_max)
		calib_error("Insufficient calibration data");
	
	/* Apply the calibration to the measured resistance */
	a_coef = calib_steinhart_hart.a_coef;
	b_coef = calib_steinhart_hart.b_coef;
	c_coef = calib_steinhart_hart.c_coef;
	R = meas_res;
	Tinv = a_coef + b_coef * log(R) + c_coef * pow(log(R), 3);
	Tk = 1. / Tinv;
	*temp = Tk - 273.15;
	
	/* Done */
	return;
}

/*	**************************************************************	*/
/*
 *	calib_plot_steinhart_hart
 *	This function plots the results of the Steinhart-Hart equation.
 *
 *	Input:
 *	calib_steinhart_hart, structure containing the Steinhart_Hart parameters
 *
 *	Zachary Wolf
 *	12/13/02
 */

void calib_plot_steinhart_hart(struct calib_steinhart_hart_struct calib_steinhart_hart)
{
	/* Declare variables */
	int num_meas = 201;
	double meas_val[201];
	double calib_out[201];
	double range_min, range_max;
	double delta_meas_val;
	int i;
	FILE* file_ptr;
	
	/* Generate a number of points within the range of the calibration data */
	range_min = calib_steinhart_hart.meas_range_min;
	range_max = calib_steinhart_hart.meas_range_max;
	delta_meas_val = (range_max - range_min) / (num_meas - 1);
	for (i = 0; i < num_meas; i++) meas_val[i] = range_min + i * delta_meas_val; 
	
	/* Apply the calibration function to the generated measured points */
	for (i = 0; i < num_meas; i++) calib_apply_calib_steinhart_hart(calib_steinhart_hart, meas_val[i], &calib_out[i]);
 
	/* Plot the results */
	/* Open the plt file */
	file_ptr = fopen("Val_steinhart_hart.plt", "w");
	if (file_ptr == NULL) return;
	
	/* Write the values to the plt file */
	fprintf(file_ptr, "%% meas_val (R ohm), calib_out (T degC)\n");
	for (i = 0; i < num_meas; i++)
		fprintf(file_ptr, "%f  %f\n", meas_val[i], calib_out[i]);
	
	/* Close the plt file */
	fclose(file_ptr);
	
	/* Done */
	return;
}

/*	**************************************************************	*/
/*
 *	calib_test_steinhart_hart
 *	This function performs tests of the Steinhart-Hart functions.
 *
 *	Zachary Wolf
 *	12/13/02
 */

void calib_test_steinhart_hart(void)
{
	/* Declare variables */
	struct calib_steinhart_hart_struct calib_steinhart_hart;
	struct calib_data_struct meas_data;
	int i;
	double calib_out[CALIB_MAX_NUM_CALIB_PTS];
	double error[CALIB_MAX_NUM_CALIB_PTS];
	FILE* file_ptr;
	
	/* Read the fake calibration data from the file */
	calib_get_calib_steinhart_hart("Val_test_steinhart_hart_file.dat", &calib_steinhart_hart);
	
	/* Check the probe name */
	if (strstr(calib_steinhart_hart.probe_name, "test_calib") == NULL) printf("Problem with the calibration data.\n");
	
	/* Plot the functional form of the calibration data */
	calib_plot_steinhart_hart(calib_steinhart_hart);
	
	/* Done */
	return;
}

/*	**************************************************************	*/

/*	GENERATE TEST CALIB DATA  */

/*	**************************************************************	*/
/*
 *	calib_generate_test_calib_data
 *	This function creates fake calibration data.
 *	This is used to test the calibration procedure.
 *
 *	Input:
 *	num_calib_pts, number of calibration points
 *	meas_range_min, minimum of the measured calibration points
 *	meas_range_max, maximun of the measured calibration points
 *
 *	Output:
 *	calib_data, structure containing the fake calibration data
 *
 *	Zachary Wolf
 *	12/13/02
 */

void calib_generate_test_calib_data(int num_calib_pts, double meas_range_min, double meas_range_max, struct calib_data_struct* calib_data)
{
	/* Declare variables */
	int i;
	double meas_calib, real_calib;

	/* Name of the probe */
	strcpy(calib_data->probe_name, "test_calib");
	
	/* Number of calibration points */
	calib_data->num_calib_pts = num_calib_pts;
	
	/* Measured data */
	for (i = 0; i < num_calib_pts; i++)
	{
		meas_calib = meas_range_min + i * (meas_range_max - meas_range_min) / (num_calib_pts - 1);
		calib_data->meas_calib[i] = meas_calib;
	}
	
	/* Real data */
	for (i = 0; i < num_calib_pts; i++)
	{
		meas_calib = calib_data->meas_calib[i];
		real_calib = meas_calib + 0.1 * pow(meas_calib, 3);
		calib_data->real_calib[i] = real_calib;
	}
	
	/* Done */
	return;
}

/*	**************************************************************	*/

/*	INTERNAL FUNCTIONS  */

/*	**************************************************************	*/
/*
 *	calib_error
 *	This function handles errors for the CALIB module.
 *
 *	Input:
 *	message, string to display in standard I/O
 *
 *	Zachary Wolf
 *	10/11/98
 */
 
void calib_error(char* message)
{
	/* Declare variables */
	char buf[CALIB_MAX_NUM_CHAR];
	
	/* Notify the operator of the error */
	printf("\nCALIB ERROR: %s\n", message);
	Beep();
	Delay(.5);
	Beep();
	
	/* Terminate the program if the operator desires */
	printf("Press ENTER to continue.\nPress any key then ENTER to terminate program.\n");
	fgets(buf, CALIB_MAX_NUM_CHAR, stdin);
	if (buf[0] == '\n') return;
	else exit(0);
}

/*	**************************************************************	*/
/*
 *	calib_message
 *	This function handles messages for the IMAG module.
 *
 *	Input:
 *	message, string to display in standard I/O
 *
 *	Zachary Wolf
 *	10/12/98
 */

void calib_message(char* message)
{
	/* Print the message */
	printf("%s\n", message);
	
	/* Done */
	return;
}

/*	**************************************************************	*/
/*
 *	calib_sort_calib_data
 *	This function sorts the calibration data for the probe so the
 *	points are ordered from the smallest measured value to the largest.
 *
 *	Input/Output:
 *	calib_data, structure containing the calibration data
 *
 *	Zachary Wolf
 *	12/13/02
 */

void calib_sort_calib_data(struct calib_data_struct* calib_data)
{
	/* Declare variables */
	int i, j;
	double min;
	int index;
	double temp;
	
	/* Sort */
	for (i = 0; i < calib_data->num_calib_pts; i++)
	{
		min = 1000.;
		for (j = i; j < calib_data->num_calib_pts; j++)
		{
			if (calib_data->meas_calib[j] < min)
			{
				min = calib_data->meas_calib[j];
				index = j;
			}
		}
		temp = calib_data->meas_calib[i];
		calib_data->meas_calib[i] = calib_data->meas_calib[index];
		calib_data->meas_calib[index] = temp;
		temp = calib_data->real_calib[i];
		calib_data->real_calib[i] = calib_data->real_calib[index];
		calib_data->real_calib[index] = temp;
	}
	
	/* Done */
	return;
}