function varargout = FitData6(varargin)
% FITDATA6 Application M-file for FitData6.fig
%    FIG = FITDATA6 launch FitData6 GUI.
%    FITDATA6('callback_name', ...) invoke the named callback.

% Last Modified by GUIDE v2.5 29-Aug-2008 08:57:56


if nargin<1  % LAUNCH GUI

    % Read Important parameters from FitParameters
    disp('Opening file FitParameters.jfit in order to read certain key fitting parameters...')
    fid = fopen('FitParameters.jfit','r');
    if fid ~= -1 % file opened successfully
        disp(['Current values read from the file ' pwd '\FitParameters.jfit :'])
        
        scrap = fgets(fid);
        workdir = fgetl(fid);
        
        scrap = fgets(fid);
        maxn = fscanf(fid,'%f');
        
        scrap = fgets(fid);
        nstartingpoints = fscanf(fid,'%f');
        
        scrap = fgets(fid);
        maxnevals = fscanf(fid,'%f');
        
        scrap = fgets(fid);
        nsim = fscanf(fid,'%f');
        
        scrap = fgets(fid);
        numbins = fscanf(fid,'%f');
        
        scrap = fgets(fid);
        lw = fscanf(fid,'%f')
        
        scrap = fgets(fid);
        fs = fscanf(fid,'%f');
        
        scrap = fgets(fid);
        gl = fscanf(fid,'%f');

        scrap = fgets(fid);
        strfittype = fscanf(fid,'%s');

        fclose(fid);
        disp('******************************************************************************')
    else 
        disp('Could not open FitParameters.jfit. Check that the file exists in the same directory as the fitting program itself.')
    end
    
    % Check values read from file and write defaults if unsuitable:
    if ~exist('workdir') || ~ischar(workdir) || ~exist(workdir,'dir')  
        msg = 'Invalid working directory read from FitParameters.jfit';
        disp(msg)
        workdir = pwd;
    end	
    
    if ~exist('maxn') || ~isnumeric(maxn) || maxn<1 || maxn>100 || isempty(maxn)
        msg = 'Invalid maximum n read from FitParameters.jfit';
        disp(msg)
        maxn=10;
    end 
    disp(['Maximum value of n which will be considered is ' num2str(maxn,'%2.0i') ]);
    disp('              Increase this if you think higher n are likely.')
    
    if ~exist('nstartingpoints') || ~isnumeric(nstartingpoints)  || nstartingpoints<1 || nstartingpoints>100 || isempty(nstartingpoints)
        msg = 'Invalid number of starting points read from FitParameters.jfit';
        disp(msg)
        nstartingpoints = 10;
    end 
    disp(['Fitting will be repeated with ' num2str(nstartingpoints,'%2.0i') ' different starting points.'])
    disp('              Increase this to make the fitting more reliable.')
    
    if ~exist('maxnevals') || ~isnumeric(maxnevals) || maxnevals<1 || maxnevals>1e20 || isempty(maxnevals)
        msg = 'Invalid number of starting points read from FitParameters.jfit';
        disp(msg)
        maxnevals = 1e6;
    end     
    options = optimset('MaxFunEvals',maxnevals);
    disp(['The program will give up fitting after ' num2str(optimget(options,'MaxFunEvals'),'%2.0f') ' function evaluations.'])
    disp('              You could try increasing this if the fitting fails a lot.')
    
    if ~exist('nsim') || ~isnumeric(nsim) || nsim>1e10 || isempty(nsim)
        msg = 'Invalid number of simulations read from FitParameters.jfit';
        disp(msg)
        nsim = 1000;
    end      
    disp(['Number of runs used in Monte Carlo simulations is ' num2str(nsim,'%5.0i') ]);
    disp('              Increase this to get more accurate distributions of measures such as chi-squared.')
    
    disp('Decreasing any of these parameters will speed up the run time, generally at the cost of accuracy.')
    disp('Display options:')
    
    if ~exist('numbins') || ~isnumeric(numbins) || numbins<1 || numbins>100 || isempty(numbins)
        msg = 'Invalid number of bins read from FitParameters.jfit';
        disp(msg)
        numbins = 30;
    end 
    disp([' Amplitude histograms will have ' num2str(numbins) ' bins.'])
    
    if ~exist('lw') || ~isnumeric(lw) || lw<0 || lw>100 || isempty(lw)
        msg = 'Invalid linewidth read from FitParameters.jfit';
        disp(msg)
        lw = 2;
    end 
    disp([' Graphs will have a linewidth of ' num2str(lw)'.'])
    
    if ~exist('fs') || ~isnumeric(fs) || fs<0 || fs>100 || isempty(fs)
        msg = 'Invalid fontsize read from FitParameters.jfit';
        disp(msg)
        fs = 12;
    end 
    disp([' Axes will have a fontsize of ' num2str(fs) '.'])
    
    if ~exist('gl') || ~isnumeric(gl) || gl<0 || gl>1 || isempty(gl)
        msg = 'Invalid glitchrate read from FitParameters.jfit';
        disp(msg)
        gl = 0;
    end 
    disp(sprintf('Fit will assume a glitch rate of %f\n',gl))
    
    if ~exist('strfittype') || ~( strcmp(strfittype,'MLE') || strcmp(strfittype,'GUY1') || strcmp(strfittype,'GUY2') ) 
        msg = 'Invalid fit-type read from FitParameters.jfit';
        disp(msg)
        strfittype = 'MLE';
    end 
    switch strfittype
        case 'MLE'            
            fittype=1;
            maximand = 'log(L)'; % for purpose of printing out info
            disp('Fit will maximise log-likelihood, ie SUM(log(p))')
        case 'GUY1'
            fittype=2;
            maximand = 'SUM(p.log(p))'; % for purpose of printing out info
            disp('Fit will maximise SUM(p.log(p))')
        case 'GUY2'
            fittype=3;
            maximand = 'SUM(p)'; % for purpose of printing out info
            disp('Fit will maximise SUM(p)')
    end

    fig = openfig(mfilename,'reuse');
    handles = guihandles(fig);
    handles.workdir = workdir;
    set(handles.datadir,'string',workdir);
    % Since we have changed the working directory, update this:
    datadir_Callback(handles.datadir, [], handles)

    % Write in these values to the fitspecs structure
    fitspecs.maxn = maxn;
    fitspecs.nstartingpoints = nstartingpoints;
    fitspecs.maxnevals = maxnevals;
    fitspecs.nsim = nsim;
    fitspecs.numbins = numbins;
    fitspecs.lw = lw;
    fitspecs.fs = fs;
    fitspecs.glitchrate = gl;
    fitspecs.fittype = fittype;
    fitspecs.maximand = maximand; % for purpose of printing out info
    
    % Default values
    fitspecs.defaultn = 4;
    fitspecs.defaultp = 0.5;
    fitspecs.defaultq = 200;
    fitspecs.defaultqs = 0;
    fitspecs.defaultqv = 1;
    fitspecs.defaultnoise = 50;
    fitspecs.defaultpcond = 1;
    fitspecs.defaultoffset = 0;
    % Resampling default:
    fitspecs.nresamp = 10;
    fitspecs.havefit = 0; % No fit is available at startup
        
    % In the absence of other info, set fit params to unconstrained
    fitspecs.n = NaN;
    fitspecs.p = NaN;
    fitspecs.q = NaN;
    fitspecs.qs = NaN;
    fitspecs.qvtype = NaN;
    fitspecs.noise = NaN;
    fitspecs.pcond = NaN;
    fitspecs.offset = NaN;
    
    handles.fitspecs = fitspecs;    
    
    if nargout > 0
        varargout{1} = fig;
    end
    
    % Write data about fit:
    set(handles.info2,'string',sprintf('Fit type = %s; glitch rate = %5.2f',strfittype,gl))
    
    set(handles.npopup,'Visible','off')
    set(handles.pedit,'Visible','off')
    set(handles.qedit,'Visible','off')
    set(handles.qsedit,'Visible','off')
    set(handles.qvpopup,'Visible','off')
    set(handles.pcondedit,'Visible','off')
    set(handles.offsetedit,'Visible','off')
    set(handles.percentsign,'Visible','off')
    
    set(handles.nbinsedit,'string',handles.fitspecs.numbins);
    set(handles.lwedit,'string',handles.fitspecs.lw);
    set(handles.fsedit,'string',handles.fitspecs.fs);
    
    set(handles.maxnedit,'string',handles.fitspecs.maxn);
    set(handles.nstartsedit,'string',handles.fitspecs.nstartingpoints);
    set(handles.maxnevalsedit,'string',handles.fitspecs.maxnevals);
    
     % Add data re the fit to the figure
    set(handles.fitpdftext,'visible','off');
    set(handles.fitcdftext,'visible','off');
   
    set(handles.qvpopup,'Visible','off','string',{'Type I';'flat (mixed Type I/II)'},'value',handles.fitspecs.defaultqv)
    
    guidata(fig, handles);



elseif ischar(varargin{1}) % INVOKE NAMED SUBFUNCTION OR CALLBACK

	try
		if (nargout)
			[varargout{1:nargout}] = feval(varargin{:}); % FEVAL switchyard
		else
			feval(varargin{:}); % FEVAL switchyard
		end
	catch
		disp(lasterr);
	end

end


%| ABOUT CALLBACKS:
%| GUIDE automatically appends subfunction prototypes to this file, and 
%| sets objects' callback properties to call them through the FEVAL 
%| switchyard above. This comment describes that mechanism.
%|
%| Each callback subfunction declaration has the following form:
%| <SUBFUNCTION_NAME>(H, EVENTDATA, HANDLES, VARARGIN)
%|
%| The subfunction name is composed using the object's Tag and the 
%| callback type separated by '_', e.g. 'slider2_Callback',
%| 'figure1_CloseRequestFcn', 'axis1_ButtondownFcn'.
%|
%| H is the callback object's handle (obtained using GCBO).
%|
%| EVENTDATA is empty, but reserved for future use.
%|
%| HANDLES is a structure containing handles of components in GUI using
%| tags as fieldnames, e.g. handles.figure1, handles.slider2. This
%| structure is created at GUI startup using GUIHANDLES and stored in
%| the figure's application data using GUIDATA. A copy of the structure
%| is passed to each callback.  You can store additional information in
%| this structure at GUI startup, and you can change the structure
%| during callbacks.  Call guidata(h, handles) after changing your
%| copy to replace the stored original so that subsequent callbacks see
%| the updates. Type "help guihandles" and "help guidata" for more
%| information.
%|
%| VARARGIN contains any extra arguments you have passed to the
%| callback. Specify the extra arguments by editing the callback
%| property in the inspector. By default, GUIDE sets the property to:
%| <MFILENAME>('<SUBFUNCTION_NAME>', gcbo, [], guidata(gcbo))
%| Add any extra arguments after the last argument, before the final
%| closing parenthesis.


% --------------------------------------------------------------------
function varargout = nconst_Callback(h, eventdata, handles, varargin)

if (get(h,'Value') == get(h,'Max'))
    % then checkbox is checked-take approriate action
    set(h,'string','constrained')
    possiblen = [1:handles.fitspecs.maxn];
    indx = find(possiblen == handles.fitspecs.defaultn);
    set(handles.npopup,'Visible','on','string',possiblen,'value',indx)
    handles.fitspecs.n = handles.fitspecs.defaultn;
    
else
    % checkbox is not checked-take approriate action
    set(h,'string','unconstrained')
    set(handles.npopup,'Visible','off')
    handles.fitspecs.n=NaN;
end

guidata(handles.figure1,handles);

% --------------------------------------------------------------------
function varargout = npopup_Callback(h, eventdata, handles, varargin)
% Set n to chosen value
handles.fitspecs.n = get(h,'value');
disp('n constrained to')
handles.fitspecs.n
guidata(handles.figure1,handles);


% --------------------------------------------------------------------
function varargout = pconst_Callback(h, eventdata, handles, varargin)
if (get(h,'Value') == get(h,'Max'))
    % then checkbox is checked-take approriate action
    set(h,'string','constrained')
    set(handles.pedit,'Visible','on','string',handles.fitspecs.defaultp)
    handles.fitspecs.p = handles.fitspecs.defaultp;
    
else
    % checkbox is not checked-take approriate action
    set(h,'string','unconstrained')
    set(handles.pedit,'Visible','off')
    handles.fitspecs.p=NaN;
end

guidata(handles.figure1,handles);

% --------------------------------------------------------------------
function varargout = pedit_Callback(h, eventdata, handles, varargin)
p = str2num(get(h,'string'))
if ~isempty(p) & isnumeric(p) & p >=0 & p <=1
    handles.fitspecs.p=str2num(get(h,'string'));
else
    errordlg('Invalid value entered for release probability.');
    set(h,'string',num2str(handles.fitspecs.p));
end

guidata(handles.figure1,handles);


% --------------------------------------------------------------------
function varargout = qconst_Callback(h, eventdata, handles, varargin)
if (get(h,'Value') == get(h,'Max'))
    % then checkbox is checked-take approriate action
    set(h,'string','constrained')
    set(handles.qedit,'Visible','on','string',handles.fitspecs.defaultq)
    handles.fitspecs.q = handles.fitspecs.defaultq;
    
else
    % checkbox is not checked-take approriate action
    set(h,'string','unconstrained')
    set(handles.qedit,'Visible','off')
    handles.fitspecs.q=NaN;
end

guidata(handles.figure1,handles);



% --------------------------------------------------------------------
function varargout = qedit_Callback(h, eventdata, handles, varargin)
q = str2num(get(h,'string'))
if ~isempty(q) & isnumeric(q) & q >=0 
    handles.fitspecs.q=str2num(get(h,'string'));
else
    errordlg('Invalid value entered for quantal size.');
    set(h,'string',num2str(handles.fitspecs.q));
end

guidata(handles.figure1,handles);


% --------------------------------------------------------------------
function varargout = qsconst_Callback(h, eventdata, handles, varargin)
if (get(h,'Value') == get(h,'Max'))
    % then checkbox is checked-take approriate action
    set(h,'string','constrained')
    set(handles.qsedit,'Visible','on','string',handles.fitspecs.defaultqs)
    handles.fitspecs.qs = handles.fitspecs.defaultqs;
    set(handles.percentsign,'Visible','on')
    
else
    % checkbox is not checked-take approriate action
    set(h,'string','unconstrained')
    set(handles.qsedit,'Visible','off')
    set(handles.percentsign,'Visible','off')
    handles.fitspecs.qs = NaN;
end

guidata(handles.figure1,handles);

% --------------------------------------------------------------------
function varargout = qsedit_Callback(h, eventdata, handles, varargin)
qs = str2num(get(h,'string'))
if ~isempty(qs) & isnumeric(qs) & qs >=0 
    handles.fitspecs.qs=str2num(get(h,'string'));
else
    errordlg('Invalid value entered for SD of quantal variance.');
    set(h,'string',num2str(handles.fitspecs.qs));
end

guidata(handles.figure1,handles);

% --------------------------------------------------------------------
function varargout = qvconst_Callback(h, eventdata, handles, varargin)


if (get(h,'Value') == get(h,'Max'))
    % then checkbox is checked-take approriate action
    set(h,'string','constrained')
    set(handles.qvpopup,'Visible','on','string',{'Type I';'flat (mixed Type I/II)'},'value',handles.fitspecs.defaultqv)
    handles.fitspecs.qvtype = handles.fitspecs.defaultqv;    
else
    % checkbox is not checked-take approriate action
    set(h,'string','unconstrained')
    set(handles.qvpopup,'Visible','off')
    handles.fitspecs.qvtype=NaN;
end

guidata(handles.figure1,handles);



% --------------------------------------------------------------------
function varargout = qvpopup_Callback(h, eventdata, handles, varargin)
handles.fitspecs.qvtype = get(h,'value');
guidata(handles.figure1,handles);

% --------------------------------------------------------------------
function varargout = noiseconst_Callback(h, eventdata, handles, varargin)

if (get(h,'Value') == get(h,'Max'))
    % then checkbox is checked-take approriate action
    set(h,'string','constrained')
    set(handles.noisedit,'Visible','on','string',handles.fitspecs.defaultnoise)
    set(handles.noisefilebutton,'Visible','on');
    handles.fitspecs.noise = handles.fitspecs.defaultnoise;
    
else
    % checkbox is not checked-take approriate action
    set(h,'string','unconstrained')
    set(handles.noisedit,'Visible','off')
    set(handles.noisefilebutton,'Visible','off');
    handles.fitspecs.pcond=NaN;
end

guidata(handles.figure1,handles);


% --------------------------------------------------------------------
%function varargout = noisedit_Callback(h, eventdata, handles, varargin)
function varargout = edit12_Callback(h, eventdata, handles, varargin)
noise = str2num(get(h,'string'))
if ~isempty(noise) & isnumeric(noise) & noise >=0 
    handles.fitspecs.noise=noise;
else
    errordlg('Invalid value entered for noise.');
    set(h,'string',num2str(handles.fitspecs.noise));
end

guidata(handles.figure1,handles);

% --------------------------------------------------------------------
function varargout = noisefilebutton_Callback(h, eventdata, handles, varargin)
% Read in a noisefile, fit a Gaussian to it
disp('noise')
handles
[handles.fitspecs.noise,handles.workdir] = fitnoise(handles);

set(handles.noisedit,'string',num2str(handles.fitspecs.noise));
guidata(handles.figure1,handles);


% --------------------------------------------------------------------
function varargout = pcondconst_Callback(h, eventdata, handles, varargin)

if (get(h,'Value') == get(h,'Max'))
    % then checkbox is checked-take approriate action
    set(h,'string','constrained')
    set(handles.pcondedit,'Visible','on','string',handles.fitspecs.defaultpcond)
    handles.fitspecs.pcond = handles.fitspecs.defaultpcond;
    
else
    % checkbox is not checked-take approriate action
    set(h,'string','unconstrained')
    set(handles.pcondedit,'Visible','off')
    handles.fitspecs.pcond=NaN;
end

guidata(handles.figure1,handles);




% --------------------------------------------------------------------
function varargout = pcondedit_Callback(h, eventdata, handles, varargin)
p = str2num(get(h,'string'))
if ~isempty(p) & isnumeric(p) & p >=0 & p <=1
    handles.fitspecs.pcond=p;
else
    errordlg('Invalid value entered for conductance probability.');
    set(h,'string',num2str(handles.fitspecs.pcond));
end

guidata(handles.figure1,handles);


% --------------------------------------------------------------------
function varargout = offsetconst_Callback(h, eventdata, handles, varargin)
if (get(h,'Value') == get(h,'Max'))
    % then checkbox is checked-take approriate action
    set(h,'string','constrained')
    set(handles.offsetedit,'Visible','on','string',handles.fitspecs.defaultoffset)
    handles.fitspecs.offset = handles.fitspecs.defaultoffset;
    
else
    % checkbox is not checked-take approriate action
    set(h,'string','unconstrained')
    set(handles.offsetedit,'Visible','off')
    handles.fitspecs.offset=NaN;
end

guidata(handles.figure1,handles);



% --------------------------------------------------------------------
function varargout = offsetedit_Callback(h, eventdata, handles, varargin)
offset = str2num(get(h,'string'));
if ~isempty(offset) & isnumeric(offset) 
    handles.fitspecs.offset=offset;
else
    errordlg('Invalid value entered for conductance probability.');
    set(h,'string',num2str(handles.fitspecs.offset));
end






% --------------------------------------------------------------------
function varargout = nbinsedit_Callback(h, eventdata, handles, varargin)

numbins = str2num(get(handles.nbinsedit,'string'));
if isempty(numbins) | numbins < 2 | numbins > 1000
    errordlg('Invalid number of bins.');
    set(handles.nbinsedit,'string',handles.fitspecs.numbins);
else
    % Redraw plot
    handles.fitspecs.numbins = numbins;
    axes(handles.pdfaxes);hold off
    guidata(handles.figure1,handles);
    PlotData(handles)
    PlotModel(handles,'r')
    % Update default prefs
    WriteFitParameters(handles)
    
end

guidata(handles.figure1,handles);



% --------------------------------------------------------------------
function varargout = lwedit_Callback(h, eventdata, handles, varargin)

lw = str2num(get(handles.lwedit,'string'));
if isempty(lw) | lw < 0 | lw > 10
    errordlg('Invalid linewidth.');
    set(handles.lwedit,'string',handles.fitspecs.lw);
else
    handles.fitspecs.lw = lw;
    hold off
    guidata(handles.figure1,handles);
    PlotData(handles);
    PlotModel(handles,'r');
    % Update default prefs
    WriteFitParameters(handles)
    
    guidata(handles.figure1,handles);

end



% --------------------------------------------------------------------
function varargout = fsedit_Callback(h, eventdata, handles, varargin)

fs = str2num(get(handles.fsedit,'string'));
if isempty(fs) | fs < 0 | fs > 50
    errordlg('Invalid fontsize.');
    set(handles.fsedit,'string',handles.fitspecs.fs);
else
    handles.fitspecs.fs = fs;
    set(handles.pdfaxes,'fontsize',fs);
    set(handles.cdfaxes,'fontsize',fs);
    % Update default prefs
    WriteFitParameters(handles)
    
end

guidata(handles.figure1,handles);

% --------------------------------------------------------------------
function varargout = nstartsedit_Callback(h, eventdata, handles, varargin)

ns = str2num(get(h,'string'));
if isempty(ns) | ns < 1 | ns > 10000
    errordlg('Invalid number of starting-points.');
    set(h,'string',handles.fitspecs.nstartingpoints);
else
    handles.fitspecs.nstartingpoints = ns;
    % Update default prefs
    WriteFitParameters(handles)
end
guidata(handles.figure1,handles);



% --------------------------------------------------------------------
function varargout = maxnedit_Callback(h, eventdata, handles, varargin)

n = str2num(get(h,'string'));
if isempty(n) | n < 1 | n > 100 
    errordlg('Invalid maximum n.');
    set(h,'string',handles.fitspecs.maxn);
else % Accept the given max n
    % There could be conflict if there is a current value of n:
    if isnan(handles.fitspecs.n)
        % If there is no n stored in memory, no problem:
        handles.fitspecs.maxn = n;
        % Update the pop-up so it shows 1 up to the new max value
        possiblen = [1:handles.fitspecs.maxn];
        % Change the default n, if necessary
        if handles.fitspecs.defaultn > n
            handles.fitspecs.defaultn = n;
        end
        indx = find(possiblen == handles.fitspecs.defaultn);
        set(handles.npopup,'string',possiblen,'value',indx)
    else
        % If there IS an n currently set, check that the proposed value does not exceed it
        if handles.fitspecs.n > n
            % If so, reject it
            errordlg('Current value of n is greater than this!');
            set(h,'string',handles.fitspecs.maxn);
        else
            % accept it, change the range of the pop-up and put it to the right value
            handles.fitspecs.maxn = n;
            % Change the default n, if necessary
            if handles.fitspecs.defaultn > n
                handles.fitspecs.defaultn = n;
            end
            % Update the pop-up so it shows 1 up to the new max value
            possiblen = [1:handles.fitspecs.maxn];
            indx = find(possiblen == handles.fitspecs.n);
            set(handles.npopup,'string',possiblen,'value',indx)
            
        end
    end
end
% Update default prefs
WriteFitParameters(handles)
guidata(handles.figure1,handles);





% --------------------------------------------------------------------
function varargout = maxnevalsedit_Callback(h, eventdata, handles, varargin)
ns = str2num(get(h,'string'));
if isempty(ns) | ns < 100 | ns > 1e20
    errordlg('Invalid maximum number of function evaluations.');
    set(h,'string',handles.fitspecs.maxnevals);
else
    handles.fitspecs.maxnevals = ns;
    % Update default prefs
    WriteFitParameters(handles)
end

guidata(handles.figure1,handles);



% --------------------------------------------------------------------
function varargout = exportbutton_Callback(h, eventdata, handles, varargin)
f=figure('PaperPositionMode','auto','name','Model/data PDF','visible','off');
a=axes;
exporthandles = handles;
exporthandles.figure1 = f;
exporthandles.pdfaxes = a;
PlotModel(exporthandles,'r');
exporthandles = guidata(exporthandles.figure1);
negloglike = NegLogLikelihood(handles.fitspecs.data,handles.model);
axes(exporthandles.pdfaxes);
title(['-ln(L) = ' num2str(negloglike)])
xlabel('EPSP amplitude')
ylabel('PDF')
% Output to jpg
jpgfile = [handles.workdir '\results.' handles.fitspecs.datafile '.fitgraph.jpg'];
print(f,'-djpeg',jpgfile);
close(f)
% Write model PDF to fit
pdffile = [handles.workdir '\results.' handles.fitspecs.datafile '.modelPDF'];
fid = fopen(pdffile,'w');
if fid == -1
    errordlg(['Unexpected error: Cannot open the file ' [resultsdir resultsfile '.modelPDF'] ' for writing results.'])
else    
    fprintf(fid,'%f %f\n', [ exporthandles.fitspecs.v ; exporthandles.fitspecs.PDF ] );
    msgbox(['Amplitude histogram graph has been exported to ' jpgfile  '. Model PDF has been written to ' pdffile '.']);
end
fclose(fid);

% --------------------------------------------------------------------
function varargout = overlaybutton_Callback(h, eventdata, handles, varargin)

if (get(h,'Value') == get(h,'Max'))
  set(h,'string','overlay fits');
elseif (get(h,'Value') == get(h,'Min'))
  set(h,'string','erase fits');
  PlotData(handles);
end




% --------------------------------------------------------------------
function varargout = fitbutton_Callback(h, eventdata, handles, varargin)
disp('FITTING==========')
% If there was already a fit, and we want overlay, redraw it in green
if handles.fitspecs.havefit == 1 & get(handles.overlaybutton,'value')==1
    hold on
    PlotModel(handles,'g');
end
set([handles.info1 handles.exportbutton handles.resampbutton handles.nresamp],'visible','off');
set(handles.waitsign,'visible','on','fontsize',16)
disp('n const to')
handles.fitspecs.n
handles.model = Fitting(handles);
set(handles.waitsign,'visible','off')
handles.fitspecs.havefit = 1;
guidata(h,handles);

PlotModel(handles,'r');
% Reveal the "Export Fit" button:
set([handles.exportbutton handles.resampbutton handles.nresamp],'visible','on')
handles = guidata(handles.figure1);

% Display model params
set(handles.text37,'visible','on','string','Fitted values');
set(handles.text38,'visible','on','string','Fitted values');
set(handles.nbest,'string',num2str(handles.model.n),'visible','on');
set(handles.pbest,'string',num2str(handles.model.p,'%5.3f'),'visible','on');
set(handles.qbest,'string',num2str(handles.model.q,'%5.3f'),'visible','on');
set(handles.qsbest,'string',num2str(handles.model.qsigma/handles.model.q*100,'%3.0f'),'visible','on');
if handles.model.QVtype == 1
    set(handles.qvbest,'string','I','visible','on');
elseif handles.model.QVtype == 2
    set(handles.qvbest,'string','flat','visible','on');
end
set(handles.noisebest,'string',num2str(handles.model.noise,'%5.3f'),'visible','on');
set(handles.pcondbest,'string',num2str(handles.model.pstim,'%5.3f'),'visible','on');
set(handles.offsetbest,'string',num2str(handles.model.offset,'%5.3f'),'visible','on');
negloglike = NegLogLikelihood(handles.fitspecs.data,handles.model);
set(handles.text23,'string',num2str(negloglike),'visible','on');
set(handles.exportbutton,'visible','on');
set(handles.testbutton,'visible','on');
guidata(handles.figure1,handles);

% --------------------------------------------------------------------
function varargout = exitbutton_Callback(h, eventdata, handles, varargin)
confirmandexit



% --------------------------------------------------------------------
function varargout = pushbutton4_Callback(h, eventdata, handles, varargin)
%function varargout = testbutton_Callback(h, eventdata, handles, varargin)
pfail = getpfail;
m = msgbox('Testing adequacy of fit... please wait','Adequacy tests','modal');

resultsfile = [handles.workdir '\results.' handles.fitspecs.datafile];
% WOrk out whether model is peaky
fid = fopen(resultsfile,'a');
npeaks = FindNPeaks(handles.model,min(handles.fitspecs.data),max(handles.fitspecs.data));
if npeaks <= 2
    fprintf(fid,'%s %d %s\n','0: model is NOT PEAKY (fitted PDF has ',npeaks,' <=2 maxima within the data range');
elseif npeaks > 2
    fprintf(fid,'%s %d %s\n','1: model is PEAKY (fitted PDF has ',npeaks,' >2 maxima within the data range');
end
fclose(fid);
handles.adeq = Adequacy2(resultsfile,handles.model,handles.fitspecs.data,pfail,handles.fitspecs.nsim);
guidata(h,handles)
try
    close(m)
catch
end
m=msgbox(['Adequacy plots have been exported to ' resultsfile '.adeqgraph1.jpg and -2.jpg.'],'Adequacy tests','modal');

function varargout = datafile_Callback(h, eventdata, handles, varargin)
% The selected data file has been changed - update the fitspecs structure
% accordingly.
jfile = get(h,'value');
files = get(h,'string');
datafile = files{jfile};

% Read in data from the new data file.
% Catch an error involving an invalid file:
fid = fopen( [handles.workdir '\' datafile],'r');
if fid == -1
    validdatafile=0;
    errordlg(['Cannot open the file ' datafile '. Try another file.'],'Datafile error','modal')
    return
else
    data = fscanf(fid,'%f');
    ndata = length(data);
    % NB some of my functions need data to be a row not a column vector. Next few lines ensure this:
    [tmp1,tmp2] = size(data);
    if tmp1>tmp2, data = data', end

    if ndata>1 
        handles.fitspecs.data = sort(data);
        handles.fitspecs.ndata = ndata;
    else
        errordlg(['The file ' datafile ' does not contain valid data. Enter another file.'],'Datafile error','modal')
        return
    end
end
handles.fitspecs.datafile = datafile;
% Remove model, since the data has changed:
if isfield(handles,'model')
    handles = rmfield(handles,'model');
end
% Hide the "Export fit" button - it will reappear when a new fit has been
% generated
set([handles.exportbutton handles.resampbutton handles.nresamp],'visible','off')
% May 16th 2008 - One peculiarity of the model with glitchrate is, I shall set the range of
% glitches to icnlude the entire span of data. This means that the model
% parameters depend on the data.
handles.fitspecs.glitchmin = min(data)-1;
handles.fitspecs.glitchmax = max(data)+1;
guidata(h,handles);
PlotData(handles)
guidata(h,handles);

function varargout = datadir_Callback(h, eventdata, handles, varargin)
% The selected data directory has been changed
workdir = get(h,'string');
% Check this is value    
if ~exist('workdir') || ~ischar(workdir) || ~exist(workdir,'dir')
    msg = sprintf('Invalid working directory. Reset to %s',pwd);
    errordlg(msg)
    workdir = pwd;
    set(h,'string',workdir)
end
handles.workdir = workdir;
% Look up what files are in the specified working directory
dir_struct = dir(workdir);
% Remove any directories; we just want files:
dir_struct=dir_struct(~[dir_struct.isdir]);
[sorted_names,sorted_index] = sortrows({dir_struct.name}');
set(handles.datafile,'String',sorted_names,'Value',1)
guidata(h,handles);


% --- Executes on button press in resampbutton.
function resampbutton_Callback(h, eventdata, handles)
% hObject    handle to resampbutton (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

disp('RESAMPLING==========')

set(handles.waitsign,'visible','on','fontsize',16)
resamphandles = handles;
% Relabel handles to avoid having original fit data overwritten:
resamphandles = rmfield(resamphandles,'nbest');

resamphandles.fitspecs.datafile = [handles.fitspecs.datafile '.tmp' ]; % details of each resampling will get written here, but there is really way too much data here; it is just for debugging
ndata = length(handles.fitspecs.data);
resampfig = figure('numbertitle','off','name','Resampling...');
resamphandles.pdfaxes = subplot(1,2,1);
resamphandles.cdfaxes = subplot(1,2,2);
for jresamp=1:handles.fitspecs.nresamp
    % Generate a new set of data by resampling
    indices = unidrnd(ndata,1,ndata);
    resamphandles.fitspecs.data = sort(handles.fitspecs.data(indices)); % this picks with replacement 
    
    PlotData(resamphandles)
    PlotModel(resamphandles,'r');
    set(resampfig,'name',sprintf('Resampling run #%d out of %d',jresamp,handles.fitspecs.nresamp));
    
    resampmodel = Fitting(resamphandles);
    resamp_n(jresamp) = resampmodel.n;
    resamp_p(jresamp) = resampmodel.p;
    resamp_q(jresamp) = resampmodel.q;
    resamp_qs(jresamp) = resampmodel.qsigma;
    resamp_qsAsPercentage(jresamp) = resampmodel.qsigma./resampmodel.q*100;
    resamp_qvtype(jresamp) = resampmodel.QVtype;
    resamp_pcond(jresamp) = resampmodel.pstim;
    resamp_noise(jresamp) = resampmodel.noise;
    resamp_offset(jresamp) = resampmodel.offset;
    resamp_glitchrate(jresamp) = resampmodel.glitchrate;
end

% fIANLLY, WRITE RESULTS TO A FILE
resampfile = [handles.workdir '\' handles.fitspecs.datafile '.RESAMP' ];
fid = fopen(resampfile,'w');
if fid==-1
    errordlg(['Unexpected error: Cannot open the file ' resampfile ' for writing resampling results.'])
    return
end
% First, write 95% confidence intervals
alpha = 2.5;
fprintf(fid,'%5.3f%% confidence intervals based on %d resampling runs:\n',100-2*alpha,handles.fitspecs.nresamp)
fprintf(fid,'            Entries show errorbars (2.5%% to 97.5%% percentiles) for each fitted parameter.\n');
fprintf(fid,'%20s | %20s | %20s | %20s | %20s | %20s | %20s | %20s | %20s | %20s | %20s \n'...
    ,'Parameters:','n','p','q in mV','qsigma in mV','qsigma as a %age','QVtype','noise','pstim','offset','glitchrate');
% Recrod the original values (fit to original data):
fprintf(fid,'%20s | %20d | %20.3f | %20.3f | %20.3f | %20.3f | %20.3f | %20.3f | %20.3f | %20.3f | %20.3f \n','Original fit:',...
    handles.model.n,handles.model.p,handles.model.q,handles.model.qsigma,handles.model.qsigma/handles.model.q*100,handles.model.QVtype,...
    handles.model.noise,handles.model.pstim,handles.model.offset,handles.model.glitchrate);
% Recrod the confidence intervals:
fprintf(fid,'%20s | ','Conf Intervals:');
vars = {'n' 'p' 'q' 'qs' 'qsAsPercentage' 'qvtype' 'noise' 'pcond' 'offset' 'glitchrate'};
for jvar=1:length(vars)
    v = vars{jvar};
    value = handles.fitspecs.(strrep(v,'AsPercentage',''));
    if isnan(value)
        % then the value was fitted: get the percentiles:
        tmp = eval(['resamp_' v]);
        pmin = prctile(tmp,alpha);
        pmax = prctile(tmp,100-alpha);
        label = sprintf('%8.3f to %8.3f',pmin,pmax);
    else
        label = sprintf('(set to %5.3f)',value);
        if jvar==4
            label = '(set as a %age)';% - you don't constrain quantal vairance to a value in mV; you can only constrain it as a percentage
        end
    end
    fprintf(fid,'%20s | ',label);
end
fprintf(fid,'\n');


% Now, record the results of individual resampling runs, for completeness
fprintf(fid,'\n\n%s\n','For completeness, here the results of the individual resampling runs:')
fprintf(fid,'%20s | %20s | %20s | %20s | %20s | %20s | %20s | %20s | %20s | %20s | %20s \n'...
    ,'Resampling run:','n','p','q in mV','qsigma in mV','qsigma as a %age','QVtype','noise','pstim','offset','glitchrate');
for jresamp=1:handles.fitspecs.nresamp
    fprintf(fid,'%20d | %20.3f | %20.3f | %20.3f | %20.3f | %20.3f | %20.3f | %20.3f | %20.3f | %20.3f | %20.3f \n',jresamp,...
        resamp_n(jresamp),resamp_p(jresamp),resamp_q(jresamp),resamp_qs(jresamp),resamp_qsAsPercentage(jresamp), ...
        resamp_qvtype(jresamp),resamp_noise(jresamp),resamp_pcond(jresamp),resamp_offset(jresamp),resamp_glitchrate(jresamp));   
end
fclose(fid)

set(handles.waitsign,'visible','off')

handles.fitspecs.havefit = 1;
guidata(h,handles);

% Reveal the "Export Fit" button:
set([handles.exportbutton handles.resampbutton],'visible','on')
handles = guidata(handles.figure1);





function nresamp_Callback(h, eventdata, handles)
% Set number of resampling runs to do.
% hObject    handle to nresamp (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)
nresamp = str2num(get(h,'string'))
if ~isempty(nresamp) & isnumeric(nresamp) & nresamp >=1
    handles.fitspecs.nresamp=nresamp;
else
    errordlg('Invalid value entered for number of resamples.');
    set(h,'string',num2str(handles.fitspecs.nresamp));
end

guidata(handles.figure1,handles);


% --- Executes during object creation, after setting all properties.
function nresamp_CreateFcn(hObject, eventdata, handles)
% hObject    handle to nresamp (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    empty - handles not created until after all CreateFcns called

% Hint: edit controls usually have a white background on Windows.
%       See ISPC and COMPUTER.
if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor'))
    set(hObject,'BackgroundColor','white');
end


