caseName = 'growthHighDimCoupled';

%% execute solver
executeCase('init');
tic
executeCase(caseName);
disp('--------HERE-----------')
toc
disp('^^^^^')

%% load result
[result pbeDefinition] = loadCase(caseName)

%% plot supersatureation
[X, Y] = extractDataFromPbeResult('bulkProperty', 1, result, pbeDefinition);
figure,hold on
plot(X, Y, 'k-'), title('supersaturation');

%% plot temperature
% [X, Y] = extractDataFromPbeResult('bulkProperty', 2, result, pbeDefinition);
% figure,hold on
% plot(X, Y, 'k-'), title('temperature');


%% plot volume over time
[X, Y] = extractDataFromPbeResult('integralProperty', 'volume', result, pbeDefinition);

figure,hold on
plot(X, Y(1,:), 'k-'), plot(X, Y(1,:)-Y(2,:), 'r.'), plot(X, Y(1,:)+Y(2,:), 'r.')
title('total volume');
legend({'total volume (left axis)', 'lower error (std)', 'upper error (std)'},...
    'Location', 'NorthWest');
xlabel('time in s')
ylabel('total volume in m^3')
% error extra y-axis
axisHandle = axes;
plot(X, Y(2,:)./Y(1,:), 'k--', 'parent', axisHandle)
set(axisHandle, 'color', 'none', 'yaxisLocation', 'right', ...
    'xTick', [])
legend('coef. of variation (right axis)', 'Location', 'NorthEast')
ylabel('coefficient of variation')

%% plot mean width over time
% figure,hold on
% [X, Y] = extractDataFromPbeResult('integralProperty', 'meanWidth', result, pbeDefinition);
% plot(X, Y(1,:), 'k-'), plot(X, Y(1,:)-Y(2,:), 'r.'), plot(X, Y(1,:)+Y(2,:), 'r.')
% title('mean width');
% xlabel('time in s')
% ylabel('total width in m')
% % error extra y-axis
% axisHandle = axes;
% plot(X, Y(2,:)./Y(1,:), 'k--', 'parent', axisHandle)
% set(axisHandle, 'color', 'none', 'yaxisLocation', 'right', ...
%     'xTick', [])
% ylabel('coefficient of variation')

%% plot surface over time
% figure,hold on
% [X, Y] = extractDataFromPbeResult('integralProperty', 'surface', result, pbeDefinition);
% plot(X, Y(1,:), 'k-'), plot(X, Y(1,:)-Y(2,:), 'r.'), plot(X, Y(1,:)+Y(2,:), 'r.')
% title('surface');
% xlabel('time in s')
% ylabel('total surface in m^2')
% % error extra y-axis
% axisHandle = axes;
% plot(X, Y(2,:)./Y(1,:), 'k--', 'parent', axisHandle)
% set(axisHandle, 'color', 'none', 'yaxisLocation', 'right', ...
%     'xTick', [])
% ylabel('coefficient of variation')

%% plot shape coordinates
figure,hold on
[X, Y] = extractDataFromPbeResult('genericintegralproperty', ...
    @(H, up)(H), result, pbeDefinition);
plot(X, Y(1:7,:))
% plot(X1, Y1(1,:)-Y1(2,:), 'r.'), plot(X1, Y1(1,:)+Y1(2,:), 'r.')
% error extra y-axis
% axisHandle = axes;hold on
% plot(X, Y(2,:)./Y(1,:), 'r--', 'parent', axisHandle)
% set(axisHandle, 'color', 'none', 'yaxisLocation', 'right', ...
%     'xTick', [])
title('shape coordinates');
legend({'1', '2', '3', '4', '5', '6', '7'})

%% plot relative shape coordinates
figure,hold on
[X, Y] = extractDataFromPbeResult('genericintegralproperty', ...
    @(H, up)(H ./ H([4 4 4 4 4 4 4],:)), result, pbeDefinition);
plot(X, Y(1:7,:)*1e-9)
% plot(X1, Y1(1,:)-Y1(2,:), 'r.'), plot(X1, Y1(1,:)+Y1(2,:), 'r.')
% error extra y-axis
% axisHandle = axes;hold on
% plot(X, Y(2,:)./Y(1,:), 'r--', 'parent', axisHandle)
% set(axisHandle, 'color', 'none', 'yaxisLocation', 'right', ...
%     'xTick', [])
title('relative shape coordinates');
legend({'1', '2', '3', '4', '5', '6', '7'})

%% volume distribution
% % volume distribution density, represented by volume bins
% binLower = 0; binUpper = 2e-13; binSpacing = 3e-15;
% % binLower = 1e-5; binUpper = 9e-5; binSpacing = 3e-6;
% boundDefinitionMatrix = [binLower:binSpacing:(binUpper - binSpacing);...
%         (binLower + binSpacing):binSpacing:binUpper];
% [T, Y] = extractDataFromPbeResult('propertyDistribution', ...
%     {'volume', 'volume', boundDefinitionMatrix}, result, pbeDefinition);
% %%
% pivotPositionVector = mean(boundDefinitionMatrix, 1);
% 
% nPivot = length(pivotPositionVector);
% figure,hold on
% for iTime = 1:5:(length(result))
%     
% %     plot(result(iTime).hcMatrix(1,1), ...
% %          result(iTime).hcMatrix(2,1), ...
% %           'm.', 'MarkerSize', 20)
%     switch mod((iTime-1)/5, 4)
%         case 0, markerString = 'kx-';
%         case 1, markerString = 'bx-';
%         case 2, markerString = 'rx-';
%         case 3, markerString = 'mx-';
%     end
%     plot(pivotPositionVector, Y(1:nPivot, iTime), markerString, ...
%         'lineWidth', 1.5)
%     % error
%     plot(pivotPositionVector, ...
%          Y(1:nPivot, iTime) + Y(nPivot+(1:nPivot), iTime), ...
%          [markerString(1) '.'])
%     plot(pivotPositionVector, ...
%          Y(1:nPivot, iTime) - Y(nPivot+(1:nPivot), iTime), ...
%          [markerString(1) '.'])
% end
% title('particle volume distr. density over volume')
% ylabel('particle volume distr. density in m^3/m^3')
% xlabel('particle volume in m^3')
% % set(gca,'xscale', 'log')
