function resultsnew=biftrefine(results,BiftParams)
% C***********************************************************************
% 
% resultsnew=biftrefine(results)
%
% This part of the BIFT series allows the user to recalculate the BIFT
% using a new range of alpha and dmax.
%
% Information is first passed to the routine as the results-structure
% which has the following components
%
%
% This program calculates the IFT of SAS data
% using a Bayesian a approach developed
% by Dr. S. Hansen, of the Danish Agricultural University.
%
% It has the advantage that the ususally required Dmax and Relaxation
% parameters are found using statistics of the data
%
% Details in Hansen, S. (2000) J.Appl.Chryst 33, 1415-1421
%
% It is a matlab revision of a fortran program provided kindly by 
% Dr. S. Hansen.
%
% C***********************************************************************

if nargin == 2,
    %if parameters are passed to the file then we assume that also a value for
    %Full, NoGui is passed along
    
    % For explanations see the "else"
    nmax=BiftParams.nmax;
    ndist=BiftParams.ndist;
    ncalcmax=BiftParams.ncalcmax;
    alphastart=BiftParams.alphastart
    alphaband=BiftParams.alphaband;
    alpharange=BiftParams.alpharange;
    dmaxrange=BiftParams.dmaxrange;
    ndiam=BiftParams.ndiam;
    nalpha=BiftParams.nalpha;
    
    
    searchparam.maxit=BiftParams.searchparam.maxit; 
    searchparam.minit=BiftParams.searchparam.minit;  
    searchparam.omegaorig=BiftParams.searchparam.omegaorig; 
    searchparam.omegareduction=BiftParams.searchparam.omegareduction;
    searchparam.omegamin=BiftParams.searchparam.omegamin;
    searchparam.dotsptol=BiftParams.searchparam.dotsptol; 
    
    verboseparam.on=BiftParams.verboseparam.on;  
    verboseparam.viewstep=BiftParams.verboseparam.viewstep; 
else
    default_bift_params
end



NoGui=0
x=results.xin;
y=results.yin;
sd=results.sigmayin;
sigy=sd;

qmin=min(x);
qmax=max(x);
mtot=max(size(x));
sd=sd.^2; 

ntot0=floor(results.evix(1)/(results.r(2)-results.r(1))) % this corresponds to the number of p(r) points 

%End of Data input

% This next section tries to come up with a good center to start the
% alpha, dmax variation from. First we estimate Rg from a guinieir
% expression, then we use this to set a fixed Dmax and let do a large alpha
% variation. The best result is then used as a starting point for full
% simplex optimization.
if NoGui==0
    biftfig=figure;
    set(gcf,'Name','Contour of Previous probabilities')
    contour(results.evix,results.eviy,exp(results.eviz))
    title('Relative Probabability vs. Dmax and \alpha')
    xlabel('Dmax')
    ylabel('\alpha')
    set(gca,'Yscale','Log')
    colorbar
    h=figstatus(['Refining Full Bayesian Analysis- choose desired Dmax and \alpha range'],biftfig);
end

alphamin=min(results.eviy)
alphamax=max(results.eviy)
nalpha=length(results.eviy)
dmax=min(results.evix)
dmaxx=max(results.evix)
ndiam=length(results.evix)
dstep=(dmaxx-dmax)/(ndiam-1);

% in a later version we may do a small scanning of the terrain
% surrounding the maximum, in order to determine which alpharange and
% dmax range would be appropriate

ncalc=inf;
while ncalc>ncalcmax 
    % entering d_range
    if ~NoGui 
        
        prompt={'Enter minimum value for Dmax (A):','Enter maximum value for Dmax (A):',....
                'Enter number of D-values:','Enter minimum value for alpha (A):','Enter maximum value for Alpha (A):','Number of Alphapoints:'};
        def={num2str(dmax),num2str(dmaxx),num2str(ndiam),num2str(alphamin),num2str(alphamax),num2str(nalpha)};
        dlgTitle='Input for D-max and alpha calculation';
        lineNo=1;
        qanswer=inputdlg(prompt,dlgTitle,lineNo,def);
        
        dmax=max(str2num(qanswer{1}),0);
        dmaxx=str2num(qanswer{2});
        ndiam=floor(max(str2num(qanswer{3}),1));
        dstep=(dmaxx-dmax)/(ndiam-1);           
        alphamin=min(str2num(qanswer{4}),str2num(qanswer{5}));
        alphamax=max(str2num(qanswer{4}),str2num(qanswer{5}));
        nalpha=floor(max(str2num(qanswer{6}),1));
        
    end
    
    dmax0=dmax;
    
    %C     check number of estimates
    
    ncalc=ndiam*nalpha;
    disp(['Total number of calculations = ',num2str(ncalc)])
    
    if ~NoGui
        if(ncalc > ncalcmax) 
            disp('That is too many calculations. It will take too long')
            disp('Please change the parameters')
        end
    else
        ncalc=0 % if NoGui then we just accept the number of calculations
    end
end % end of alpha,dmax span definition
resultsnew=biftcalc(results.name,x,y,sd,sigy,NoGui,biftfig,h,nmax,alphamin,alphamax,nalpha,ndiam,dmax0,dstep,ntot0,mtot,searchparam,verboseparam);
delete(biftfig)

