image thumbnail

Multivariate Data Analysis and Monitoring for the Process Industries

by

 

26 Feb 2008 (Updated )

Files from the webinar Multivariate Data Analysis and Monitoring for the Process Industries

getCurrentProcessStateAsImage()
function imagedata = getCurrentProcessStateAsImage()

data = load('pcamodel');

numvars = size(data.varnames,2);

%% Create main figure
hfig = figure('units','normalized','position',...
    [0.1,0.1,0.8,0.6],'Visible','off');

%% Create 12 axes in the figure
for i = 1:12
    h(i) = axes('parent', hfig);
end

set(h,'units','normalized',...
    'XTick', [], 'YTick', [], 'Color', 'none',...
    'XColor', [0.5,0.5,0.5], 'YColor', [0.5,0.5,0.5],...
    'NextPlot', 'add', 'Box', 'on');

set(h,{'position'},...
    {[0.05,0.55,0.2,0.4];...
    [0.3,0.55,0.2,0.4];...
    [0.05,0.3,0.45,0.2];...
    [0.05,0.05,0.45,0.2];...
    [0.55,0.875,0.4,0.075];...
    [0.55,0.757,0.4,0.075];...
    [0.55,0.639,0.4,0.075];...
    [0.55,0.521,0.4,0.075];...
    [0.55,0.403,0.4,0.075];...
    [0.55,0.285,0.4,0.075];...
    [0.55,0.167,0.4,0.075];...
    [0.55,0.05,0.4,0.075]});

%% Add loadings to biplot
axes(h(1));
xlabel('Amount of Impurity','Color','k'); ylabel('Pattern of Impurity','Color','k'); xlim([-18,4]); ylim([-8,5]);
scatter(data.loadings(:,1)*data.scale(1),data.loadings(:,2)*data.scale(2),'ro','filled');
for i=1:numvars
    line([0,data.loadings(i,1)*data.scale(1)],...
        [0,data.loadings(i,2)*data.scale(2)],'Color','r');
    text(data.loadings(i,1)*data.scale(1)-2,data.loadings(i,2)*data.scale(2),data.varnames{i});
end

%% Add control limits and labels to polar plot
axes(h(2))
xlim([-10,10]); ylim([-10,10]);
upperlimits = [4,3,-2,-3];
lowerlimits = [3,2,-3,-4];
colors = {'r','y','y','r'};
offset = 4;
for i = 1:4
    [X,Y]=pol2cart(...
        [0:1/99:1,1:-1/99:0]*2*pi,...
        [ones(1,100)* (upperlimits(i) + offset),...
        ones(1,100)*(lowerlimits(i) + offset)]);
    patch(X,Y,colors{i},'EdgeColor','none','FaceAlpha',0.4);
end
[X,Y]=pol2cart(...
    (0:1/numvars:1)*2*pi,...
    ones(1,numvars+1)*(4 + offset));
for i=1:numvars
    htext=text(X(i),Y(i),data.varnames{i});
    circlepos = (i-1)/numvars;
    if circlepos < 0.25
        set(htext,'VerticalAlignment','bottom','HorizontalAlignment','left');
    elseif circlepos == 0.25
        set(htext,'VerticalAlignment','bottom','HorizontalAlignment','center');
    elseif circlepos < 0.5
        set(htext,'VerticalAlignment','bottom','HorizontalAlignment','right');
    elseif circlepos == 0.5
        set(htext,'VerticalAlignment','middle','HorizontalAlignment','right');
    elseif circlepos < 0.75
        set(htext,'VerticalAlignment','top','HorizontalAlignment','right');
    elseif circlepos == 0.75
        set(htext,'VerticalAlignment','top','HorizontalAlignment','center');
    elseif circlepos < 1
        set(htext,'VerticalAlignment','top','HorizontalAlignment','left');
    end
end


%% Add control limits to T2 and Residual plots

for i = [3,4]
    axes(h(i))
    xlim([0,730]);
    switch i
        case 3
            ylabel('T^{2}','Color','k');
            ylim([-5,100]);
            lims = data.tcrit;
        case 4
            xlabel('Observation Number','Color','k');
            ylabel('Residual','Color','k');
            ylim([0,6]);
            lims = data.distcrit;
    end
    upperlimits = [lims(1),lims(2)];
    lowerlimits = [lims(2),lims(3)];
    colors = {'r','y'};

    for j=1:2
        patch([0,0,730,730],...
            [upperlimits(j),lowerlimits(j),lowerlimits(j),upperlimits(j)],...
            colors{j},'EdgeColor','none','FaceAlpha',0.4);
    end
    
end


%% Add control lmits to variable plots
for i = 5:12
    axes(h(i));
    xlim([0,730]); ylim([-4,4]);
    if i==12
        xlabel('Observation Number','Color','k');
    end
    ylabel(data.varnames{i-4},'Color','k');
    upperlimits = [4,3,-2,-3];
    lowerlimits = [3,2,-3,-4];
    colors = {'r','y','y','r'};
    for j=1:4
        patch([0,0,730,730],...
            [upperlimits(j),lowerlimits(j),lowerlimits(j),upperlimits(j)],...
            colors{j},'EdgeColor','none','FaceAlpha',0.4);
    end
    
end

%% Make a connection to the database and clear the table
s.DataReturnFormat = 'numeric';
s.ErrorHandling = 'store';
s.NullNumberRead = 'NaN';
s.NullNumberWrite = 'NaN';
s.NullStringRead = 'null';
s.NullStringWrite = 'null';
s.JDBCDataSourceFile = '';
s.UseRegistryForSources = 'yes';
s.TempDirForRegistryOutput = '';
setdbprefs(s)

conn = database('MSPM','root','admin');

nodata=1;
while nodata
    curs2 = exec(conn,'SELECT ALL Ag,Ar,Bi,Ni,Pb,Sb,Se,Te FROM metals');
    curs2 = fetch(curs2);
    metals = curs2.Data;
    if isnumeric(metals)
        nodata = 0;
    end
end
close(curs2);

close(conn);

% Standardise data, and calculate PC scores, T2, and residuals
[nX,nY]=size(metals);
metalstd = (metals-repmat(data.means,nX,1))./repmat(data.sds,nX,1);
scores = metalstd*data.loadings(:,1:2);
t2 = sum((scores(:,1:2).^2)./repmat(data.variances(1:2)',nX,1),2);
residuals = (metalstd - scores*data.loadings(:,1:2)');
distance = sqrt(sum(residuals.^2,2)./(numvars-2));

% Plot PC Scores on biplot
colors = [0,linspace(0.5,0.8,27);0,linspace(0.5,0.8,27);ones(1,28)]';
if size(scores,1)<28 %TODO Left off here.
    scatter(h(1),scores(:,1),scores(:,2),50,colors(1:nX,:),'filled');
else
    scatter(h(1),scores(end-27:end,1),scores(end-27:end,2),10,colors,'filled');
end

% Plot variables on polar plot
offset = 4;
[X,Y]=pol2cart(...
   (0:1/numvars:1)*2*pi,...
   [metalstd(end,:),metalstd(end,1)] + offset);
plot(h(2),X,Y);

% Plot T2
plot(h(3),1:nX,t2);

% Plot Residual
plot(h(4),1:nX,distance);

% Plot variables

for i = 5:12
    plot(h(i),1:nX,metalstd(:,i-4));
    set(h(i),'XLim',[0 730],'YLim',[-4 4]);
end

% Make correct settings for the pixel to dpi transformation
set(hfig,'PaperUnits','inches');
resolution = get(0, 'screenp');

% redraw the image using dpi resolution
pos = get(hfig, 'position');
set(hfig,'PaperPosition', pos ./ (resolution));

% Create the output filename
% imagedata = hardcopy(hfig, '-dopengl', '-r96');
imtmp = getframe(hfig);
imagedata = imtmp.cdata;

delete(hfig);

Contact us