% spec_gen
%	This script takes the output from the slit_fn script and uses it
%	to generate a spectrometer spectrum from an input spectrum. The
%	input spectrum is binned as a function of wavelength, so we must 
%	convert it to pixels.
%
%	   * ONLY RUN THIS SCRIPT AFTER RUNNING SLIT_FN FIRST! *

%	INPUT FILES:
%		file containing the raw spectrum
%
%	VARIABLES:
%	Input variables:
%		s_fn			matrix containing slit functions, ie the output of slit_fn
%		centre		tells the centre coordinate of each slit function in s_fn
%		n_pix			# of pixels on our CCD
%		spec_in		array containing our input spectrum
%		wv_in			array containing input wavelengths (ie x-axis for spec)
%							(units = nanometers)
%		a, b, c		params for pixel/wavelength conversion
%							pixel # = p, wavelength = L
%							p = a + b*L + c*L^2
%	
%	Output and other variables:
%		spec			array containing the rebinned input spectrum	
%
%	OUTPUT:
%		a simulated spectrometer spectrum, ie the input spectrum convolved with
%		our slit function, s_fn.
%
%	COMMENTS:
%		Takes about 10 seconds to run with n_pix = 2000. (June 28)

% Initialize. Use values of a & b from A_B.REF (for lack of any better ideas):
dir_input = 'c:\james\analyse\input_files\';
fid = fopen(strcat(dir_input,'a_b.ref'),'r');
if fid == -1 	% file was not opened
   disp('SPEC_GEN.M: Error opening A_B.REF - file not opened.')
   return
end
fgetl(fid);
a = fscanf(fid,'%f20',1);
b = fscanf(fid,'%f20',1);    % a,b from absolute reference
fclose(fid);

c = 0;	% only linear relation between pixels and wavelengths for now.

j = 1:n_pix;
% Read in spectrum:
% uh...do something here. Read in values for spec_in and wv_in.
%temp:
wv_in(j) = 450;
spec_in(j) = 10;

% Convert wavelengths to pixel numbers:
pix(j) = a + b*wv_in(j) + c*wv_in(j).^2;
% BUT this probably doesn't give a 1 to 1 correspondence of wavelength to 
% pixel number right off the bat. So we need to re-bin the spectrum so as
% to get a spectrum that is a continuous function of pixel number:
spec = interp1(pix, spec_in, [1:n_pix], 'linear');

spec_out = zeros(1,n_pix);
% Here is where we need the output from SLIT_FN: s_fn, centre, and max_len
for m = 1:1:n_pix
   for n = 1:1:max_len
      if (m + centre - n >= 1) & (m + centre - n <= n_pix)
      	spec_out(m) = spec_out(m + centre - n) + spec(m).*s_fn(m, n);
   	end      
   end
end



