function show_region_areal_growth(varargin)
% function show_region_areal_growth(varargin)
%
% Shows the areal growth real and anisotropy for ungrowth data
%
% example 1
% show_region_areal_growth;
%
% example 2
% show_region_areal_growth('ungrowthfilename', 'Staged\Data\UngrowthData\Ungrowth_aih_001.mat', 'scale', 0);
%
% Dr. A. I. Hanna (2007)


if mod(length(varargin),2) ~= 0
    % input args have not com in pairs, woe is me
    error('Arguments to show_region_areal_growth must come param/value in pairs.')
end
ungrowthfilename = [];
scale = [];
for i=1:2:length(varargin)
    switch lower(varargin{i})
        case 'ungrowthfilename'
            ungrowthfilename = varargin{i+1};
        case 'scale'
            scale = varargin{i+1};
        otherwise
            error(['Unknown parameter name passed to show_region_areal_growth.  Name was ' varargin{i}])
    end
end

if isempty(ungrowthfilename)
    [filename, pathname] = uigetfile('*.mat', 'Pick an Ungrowth File');
    if isequal(filename,0) || isequal(pathname,0)
        disp('User pressed cancel')
        return;
    end
    ungrowthfilename = [pathname, filesep, filename];
end
if isempty(scale)
    scalebtn = questdlg('Scale output?', 'Scale Question', 'Yes', 'No', 'No');
    if isempty(scalebtn)
        return;
    end
    switch scalebtn,
        case 'Yes',
            scale = 1;
        case 'No',
            scale = 0;
    end % switch
end
d = load(ungrowthfilename);
if ~isfield(d, 'ungrown_grids')
    fprintf('Invalid ungrowth file\n');
    return;
end
ungrown_grids = d.ungrown_grids;
if ~isfield(ungrown_grids, 'stage_delta_t')
    fprintf('You have not specified an inter stage interval time, defaulting to 20\n');
    T = 20;
else
    T = d.stage_delta_t;
  	fprintf('Using %3.3f as the time in hours between stages\n', T);
end
gridname = d.grid_name;
ug = d.ungrown_grids;
%clear d;
%gp = gp.ungrowth_tensors;

N = length(ug);
Acols = cell(N, 1);
Gcols = cell(N, 1);
for i=2:N
    [doubling_time, anisotropy] = growthParameters2GrowthRateAnisotropy(ug{i}.ungrowth_tensors, T);
    Acols{i} = getAnisotropyColours(anisotropy);
    Gcols{i}  = getDoublingTimeColours(doubling_time);
end
Acols{1} = zeros(size(ug{1}.ungrowth_tensors,3), 3);
Gcols{1} = zeros(size(ug{1}.ungrowth_tensors,3), 3);
drawRegions(ug, Gcols, scale, 'Growth Rate');
drawRegions(ug, Acols, scale, 'Anisotropy');
return;
%%%
%
%%%
function drawRegions(ug, cols, scale, str)
figcol = [0 0 0];

M = 3;
N = 3;
indx = 1;
window_id = 1;
figh = figure;
clf(figh); set(figh, 'Color', figcol);
set(figh, 'Name', [str, ': ', num2str(window_id)]);

for i=1:length(ug)
    g = ug{i};
    if indx>9
        window_id = window_id + 1;
        figh = figure;
        set(figh, 'Name', [str, ': ', num2str(window_id)]);
        set(figh, 'Color', figcol);
        indx = 1;
    end
    subplot(M, N, indx); indx = indx + 1;
    hold on;
  
    pts = g.tri_info.pts;
    T = g.tri_info.T;
    
    r_ind = g.region_tri_ind;
    centroids = calcStageGridCentroids(g);
    num_regions = length(g.region_tri_ind);
    plot_region_patches(g.region_tri_ind, T, pts, gca, cols{i});
    for j=1:num_regions
        %x = pts(T(g.region_tri_ind{j},:),1);
        %y = pts(T(g.region_tri_ind{j},:),2);
%        if ~isempty(x)

       %     triplot(T(g.region_tri_ind{j},:), pts(:,1), pts(:,2), cols{i}, 'Parent', gca)
           mu = centroids(j,:);
       %     x = [x; x(1)];
       %     y = [y; y(1)];
       %     fill(x, y, cols{i}(j,:));
        %    plot(gca, x, y, 'Color', [0 .8 0], 'LineWidth', 1);
            if i>1
                T = g.ungrowth_tensors(:,:,j);
               % T = gtlib_growthParams2Tensor(gp{i-1}.s_max(j), gp{i-1}.s_min(j), gp{i-1}.theta(j));
                %  r = sqrt(polyarea(x, y))/5;
                r =10;
                gtlib_plotGrowthTensor('growth_tensor', T, 'offset', mu, 'colour', [0 0 0], 'Parent', gca, 'scale', r);
            end
      %  end
    end
    if scale
        if (indx==2) && (window_id==1)
            axis image ij;
            xlim = get(gca, 'XLim');
            ylim = get(gca, 'YLim');
        else
            axis image ij;
            set(gca, 'XLim', xlim, 'YLim', ylim);
        end
    else
        axis image ij;
    end
    set(gca, 'visible', 'off');
end
return;

