0;

widthCM = 13.5;
fig = figure('units','centimeters','position',...
    [5 -5 widthCM widthCM]);

dim = 2;

if dim == 2
    X_1 = 1.2;
    lambda = 2*pi/3;
    machvec = 0:0.01:0.1;
    lengthvec = (2.5:0.5:10)'/lambda;
    n_contours = 15;
else
    machvec = 0:0.0025:0.025;
    lengthvec = (0.5:0.5:2)';
    n_contours = 20;
end
machmach = repmat(machvec,length(lengthvec),1);
lengthlength = repmat(lengthvec,1,length(machvec));
s_hats_hat = zeros(length(lengthvec),length(machvec));
    

for I_L = 1:length(lengthvec)
    if dim == 2
        L = lengthvec(I_L)*lambda;
    end
    if L < 10
        lengthstring = strcat('0',sprintf('%02.1f',L));
    else
        lengthstring = sprintf('%02.1f',L);
    end
    if dim == 2
        filename = strcat('s_hat_machdata2D_',lengthstring,'_omega3.0');
        load(filename,'s_hat2vec','geomparams')
        X_i = geomparams.X_i;
        s_hats_hat(I_L,:) = (abs(s_hat2vec'/X_i));
    else
        filename = strcat('s_hat_machdata3D_',lengthstring,'_omega0.1');
        load(filename,'s_hat2vec','geomparams')
        R_i = geomparams.R_i;
        s_hats_hat(I_L,:) = (abs(s_hat2vec'/R_i));
    end
end

upperclim = max(s_hats_hat,[],'all');
lowerclim = min(s_hats_hat,[],'all');
midpoint = mean(s_hats_hat(:,1));
maxdeviation = max([upperclim-midpoint midpoint-lowerclim]);


meandrift = mean(s_hats_hat(:,end) - s_hats_hat(:,1));

if dim == 2
    contourlevels = (midpoint-maxdeviation)...
        :maxdeviation/n_contours:midpoint+maxdeviation;
elseif dim == 3
    lowerclim = min(s_hats_hat,[],'all');
    contourlevels = lowerclim:2*midpoint/n_contours:2*midpoint;
end



set_custom_colormap(1)
[~,c] = contourf(machmach,lengthlength,s_hats_hat);
c.LevelList = contourlevels;

% p = pcolor(machmach,lengthlength,s_hats_hat);
% p.EdgeColor = 'none';

clim([contourlevels(1),contourlevels(length(contourlevels))]);

xlabel('$$M$$')
ylabel('$$L/\lambda$$')
cbar = colorbar;
cbar.Label.Interpreter = 'latex';
cbar.Label.String = '$$\tilde{s}_0^1/X_1$$';
exportgraphics(gcf,'s_hat_machplot.pdf','ContentType','vector')

