%% Plot data presented in manuscript
clear
close all

subCaseDir = 'deformIndvCases';
% Load in data sets:
cMapsAe = load('contourMapsAe.mat'); % Data for Fig. 2
deformAe = load('deformAe.mat'); % Data for Fig. 3a,b
deformScale = load('deformScaling.mat'); % Data for Fig. 3c
cMapsAlpha = load('contourMapsAlpha.mat'); % Data for Fig. 4
indvCase = load('ms003mpt010.mat');


%% CL vs Ae and alpha contour plot
CL_range = [0 2.5]; % Colorbar range
nCont = 20;
figure; tiledlayout(1,2);

% Flexible Membrane Wings:
nexttile;
[~,h] = contourf(cMapsAe.AeGrid, cMapsAe.alphaGrid, cMapsAe.CLGrid, linspace(CL_range(1),CL_range(2),nCont));
clim(CL_range)
set(h,'linestyle','none'); set(gca,'XScale', 'log')
xlabel('$$Ae$$', 'interpreter', 'latex');
ylabel('$$\alpha^*$$', 'interpreter', 'latex');
title('Flexible Membrane wings')

% Rigid Wings:
nexttile;
[~,h] = contourf(cMapsAe.AeRigid, cMapsAe.alphaRigid, cMapsAe.CLGridRigid, linspace(CL_range(1),CL_range(2),nCont));
clim(CL_range)
set(h,'linestyle','none'); set(gca,'XScale', 'log')
xlabel('$$Ae$$', 'interpreter', 'latex');
ylabel('$$\alpha^*$$', 'interpreter', 'latex');
cb = colorbar; ylabel(cb, '$$\bar{C}_L$$', 'interpreter', 'latex')


%% Efficiency eta vs Ae and alpha contour plot
eta_range = [0 1]; % Colorbar range
nCont = 20;
figure; tiledlayout(1,2);

% Flexible Membrane Wings:
nexttile;
[~,h] = contourf(cMapsAe.AeGrid, cMapsAe.alphaGrid, cMapsAe.etaGrid, linspace(eta_range(1),eta_range(2),nCont));
clim(eta_range);
set(h,'linestyle','none'); set(gca,'XScale', 'log')
xlabel('$$Ae$$', 'interpreter', 'latex');
ylabel('$$\alpha^*$$', 'interpreter', 'latex');
title('Flexible Membrane wings')

% Rigid Wings:
nexttile;
[C,h] = contourf(cMapsAe.AeRigid, cMapsAe.alphaRigid, cMapsAe.etaGridRigid, linspace(eta_range(1),eta_range(2),nCont));
clim(eta_range)
set(h,'linestyle','none'); set(gca,'XScale', 'log')
xlabel('$$Ae$$', 'interpreter', 'latex');
ylabel('$$\alpha^*$$', 'interpreter', 'latex');
title('Rigid wings')
cb = colorbar; ylabel(cb, '$$\bar{\eta}$$', 'interpreter', 'latex')


%% Membrane camber VS Aeroelastic number
figure, 
semilogx(deformAe.Ae, deformAe.maxCamber, '.', 'MarkerSize', 20)
xlabel('$$Ae$$', 'interpreter', 'latex');
ylabel('$$\hat{z}_{max} / c$$', 'interpreter', 'latex');


%% Leading and trailing edge rotation angle VS Aeroelastic number
figure,
semilogx(deformAe.Ae, deformAe.maxGammaLE, '.', 'MarkerSize', 20)
hold on;
semilogx(deformAe.Ae, deformAe.maxGammaTE, '.', 'MarkerSize', 20)
xlabel('$$Ae$$', 'interpreter', 'latex');
ylabel('$$\hat{\gamma}$$', 'interpreter', 'latex');
legend('LE', 'TE')


%% Membrane camber VS Normal force coefficient scaling
figure, 
semilogx(deformScale.avgCN./(deformScale.Ae.*sind(deformScale.alpha)), ...
    deformScale.maxCamberProj, '.', 'MarkerSize', 20)
xlabel('$$\frac{We_{CN}}{\sin(\alpha)}$$', 'interpreter', 'latex');
ylabel('$$\hat{z}_{max} / c$$', 'interpreter', 'latex');


%% CL vs alpha_LE and alpha_TE contour plot
CL_range = [0 2.5]; % Colorbar range
nCont = 20;

figure;
[~,h] = contourf(cMapsAlpha.alphaTE, cMapsAlpha.alphaLE, cMapsAlpha.CLGrid, linspace(CL_range(1),CL_range(2),nCont));
set(h,'linestyle','none'); clim(CL_range); axis equal;
xlabel('$$Ae$$', 'interpreter', 'latex');
ylabel('$$\alpha^*$$', 'interpreter', 'latex');
cb = colorbar; ylabel(cb, '$$\bar{C}_L$$', 'interpreter', 'latex')


%% Efficiency vs alpha_LE and alpha_TE contour plot
eta_range = [0 1]; % Colorbar range
nCont = 20;

figure;
[~,h] = contourf(cMapsAlpha.alphaTE, cMapsAlpha.alphaLE, cMapsAlpha.etaGrid, linspace(eta_range(1),eta_range(2),nCont));
set(h,'linestyle','none'); clim(eta_range); axis equal;
xlabel('$$Ae$$', 'interpreter', 'latex');
ylabel('$$\alpha^*$$', 'interpreter', 'latex');
cb = colorbar; ylabel(cb, '$$\bar{\eta}$$', 'interpreter', 'latex')


%% Plot one individual case for the deformation and force measurements
figure(); tiledlayout(1,2);
for framei = 1:length(indvCase.tT)
    nexttile(1); hold off;
    plot(indvCase.xWing(:,framei) / indvCase.c, indvCase.yWing(:,framei) / indvCase.c, "LineWidth", 1.0, "Color", [0.53 0.80 0.53])
    hold on; axis equal;
    plot(indvCase.xWing(1:2,framei) / indvCase.c, indvCase.yWing(1:2,framei) / indvCase.c, "LineWidth", 2.0, "Color", [0.2 0.2 0.2])
    plot(indvCase.xWing(end-1:end,framei) / indvCase.c, indvCase.yWing(end-1:end,framei) / indvCase.c, "LineWidth", 2.0, "Color", [0.2 0.2 0.2])
    xlim([-0.1 0.9]); ylim([-1, 0.1])
    xlabel('$$x / c$$', 'interpreter', 'latex');
    ylabel('$$y / c$$', 'interpreter', 'latex');

    nexttile(2); hold off;
    plot(indvCase.tT, indvCase.CL, "LineWidth", 2.0)
    hold on; xlim([0 0.5])
    plot(indvCase.tT(framei), indvCase.CL(framei), 'o', "LineWidth", 2.0)
    xlabel('$$t / T$$', 'interpreter', 'latex');
    ylabel('$$C_L$$', 'interpreter', 'latex');
    title(sprintf('Ae=%.2f, AOA=%.0fdeg', indvCase.Ae, indvCase.alphaA))

    pause(0.05)
end
