https://github.com/Data2Dynamics/d2d
Raw File
Tip revision: d72d4c223fb8adc0e3e44a1d3d1eafce18b995ae authored by Andreas Raue on 06 March 2015, 00:17:07 UTC
Increase C version code, after Joep's submits today to prevent compatibility problems
Tip revision: d72d4c2
arPlotMCMCParVsChi.m
% plot mcmc chains

function arPlotMCMCParVsChi(jks, Nthinning, popt, chi2opt)

global ar

if(~exist('Nthinning','var'))
    Nthinning = 1;
end
ps = ar.ps;
chi2s = ar.chi2s;
if(Nthinning>1)
    ps = ps(mod(1:size(ps,1),Nthinning)==1,:);
    chi2s = chi2s(mod(1:length(chi2s),Nthinning)==1,:);
end

h = myRaiseFigure('mcmc chains');
set(h, 'Color', [1 1 1]);

if(~exist('jks','var') || isempty(jks))
    jks = find(ar.qFit==1);
end


jks = jks(ar.qFit(jks)==1);

[nrows, ncols] = arNtoColsAndRows(length(jks));

count = 1;
for jk=jks
    
    xlimtmp2 = (max(ps(:,jk))-min(ps(:,jk)))*0.05;
    if(xlimtmp2>0)
        xlimtmp = [min(ps(:,jk))-xlimtmp2 max(ps(:,jk))+xlimtmp2];
    end
    
    g = subplot(nrows, ncols, count);
    arSubplotStyle(g)

    plot(ps(:,jk), chi2s, 'kx', 'MarkerSize', 1)
    if(nargin>2)
        hold on
        % colors = jet(length(popt));
        for j=1:length(popt)
            % plot(popt{j}(jk), chi2opt(j), '*', 'Color', colors(j,:));
            plot(popt{j}(jk), chi2opt(j), '*r');
        end
        hold off
    end
    
    xlim([xlimtmp(1)-xlimtmp2*0.05 xlimtmp(2)+xlimtmp2*0.05]);
    
    if(nargin>2)
        ylim([min([chi2s(:); chi2opt(:)]) quantile(chi2s, 0.99)]);
    else
        ylim([min(chi2s(:)) quantile(chi2s, 0.99)]);
    end
    title(arNameTrafo(ar.pLabel{jk}))
    
    count = count + 1;
end



function h = myRaiseFigure(figname)
global pleGlobals
openfigs = get(0,'Children');

figcolor = [1 1 1];

if(isfield(pleGlobals, 'fighandel_multi') && ~isempty(pleGlobals.fighandel_multi) && ...
    pleGlobals.fighandel_multi ~= 0 && ...
    sum(pleGlobals.fighandel_multi==openfigs)>0 && ...
    strcmp(get(pleGlobals.fighandel_multi, 'Name'), figname))

    h = pleGlobals.fighandel_multi;
    figure(h);
else
    h = figure('Name', figname, 'NumberTitle','off', ...
        'Units', 'normalized', 'Position', ...
        [0.1 0.1 0.6 0.8]);
    set(h,'Color', figcolor);
    pleGlobals.fighandel_multi = h;
end



back to top