function B = extract_shape_PC_data(modeldirectory)
% function B = extract_shape_PC_data(modeldirectory)
%
% This is a simple script that takes ALL the point models related to the
% current template and stores the points, the aligned points, and the pc
% values to a matlab binary file as defined by the user.
%
% related: shapematrix2excell, exportModelData, AAMToolbox
%
% Dr. A. I. Hanna (2006).

% Get the shape model
if nargin<1
    modeldirectory = uigetdir(pwd, 'Pick a Model Directory');
    if modeldirectory == 0
        fprintf('User pressed cancel');
        return;
    end
end
[shpfilename, shppathname] = uiputfile('*.mat', 'Save As...');
if isequal(shpfilename,0) || isequal(shppathname,0)
    disp('User pressed cancel')
    return;
end
pdmfilename = [modeldirectory, filesep, 'mod_pdm.mat'];
detailsfilename = [modeldirectory, filesep, 'details.mat'];
modelfilename = [modeldirectory, filesep, 'model.aam_dat'];
aam = load(modelfilename, '-mat');
aam = aam.aam;
template = get(aam, 'PointModelTemplate');
[templatepath, templatename, ext, vers] = fileparts(get(template, 'name'));
pmdir = ['PointModels', filesep, templatename];
pdm = load(pdmfilename);
pdm = pdm.pdm;
% details = load(detailsfilename);
% details = details.details;
details = get(aam, 'procrustesParams');
% Set the PM directory
% For each group, take the first name in that group and estimate that shape
% by the model.
var_covered = pdm.FractionPCs;
Xm = pdm.Xm;
P = pdm.pca.P;
b = pdm.pca.b;
for i=1:length(b); str{i} = num2str(i); end;
[s,v] = listdlg('PromptString','Select PCs:', 'SelectionMode','multiple', 'ListString',str, 'ListSize', [300 300]);
if v == 0
    fprintf('User pressed cancel\n');
    return;
end
total_var = sum(b);
P = P(:,s);
b = b(s);
PCs = s;
dim =2;
fprintf('Variance covered: %f%%\n', 100*(sum(b)/total_var));
opts = [details.scaling, details.rotation, details.translation];
pmfiles = dir([pmdir, filesep, '*_pm.mat']);
pmfiles = {pmfiles.name};

[s,v] = listdlg('PromptString','Select PM files to extract:', 'SelectionMode','multiple', 'ListString',pmfiles, 'ListSize', [400 300]);
if isempty(s)
    fprintf('User hit cancel\n');
    return;
end
pmfiles = {pmfiles{s}};
shape_pc_values = [];
M = length(pmfiles);
waitbar_h = waitbar(0, 'Exporting shape PCs...');
for i=1:M
    waitbar(i/M, waitbar_h);
    pmfilename = [pmdir, filesep, pmfiles{i}];
    if exist(pmfilename)
        pts = load(pmfilename);
        pts = pts.pts;
        [aligned_pts, template] = pcalib_GPA('data', pts, 'opts', details, 'dimension', dim, 'template', Xm);
        model_b = P'*(aligned_pts - Xm);
        model_b_std = model_b./sqrt(b);
    else
        model_b = -1*ones(size(b));
        model_b_std = -1*ones(size(b));
    end
    shape_pc_values(i).filename = pmfilename;
    shape_pc_values(i).pc_values = model_b;
    shape_pc_values(i).std_pc_values = model_b_std;
    shape_pc_values(i).pts = pts;
    shape_pc_values(i).aligned_pts = aligned_pts;
end
save([shppathname, filesep, shpfilename], 'shape_pc_values', 'PCs');
close(waitbar_h);