% Script to plot variance in u, v, rho.
% 
% This script assumes that all relevant variables are already present in the workspace.

for i = 0 : round(Nt/app.PlotIntrvl)
%for i = 0:2

   dataname  = sprintf('%s/%s/%03d%s', d0, app.savemat, i, suffix_string);
   imagename = sprintf('%s/%s/var%03d%s.png', d0, app.saveimg, i, suffix_string);

   if ( ((exist(dataname,'file')==2) | (exist([dataname,'.mat'],'file')==2)) & ...
        ( (exist(imagename,'file')~=2) | redoPlots) )
    
      load (dataname);

      % Extract variance of u, v, rho.

%      Cyy = (Yi'*Yi);
%      [D] = eig(Cyy);extrem(D)

%      var_u = diag(ui*Cyy*ui')*((nu_x/h)^2);
%      var_v = diag(vi*Cyy*vi')*((nu_x/h)^2);
%      var_rho = diag(rhoi*Cyy*rhoi')*(DelS*DelS);
      % Create sxs covariance matrix of coefficients.
      CYY = cov(Yi);

      % Perform its SVD.
      [V,D] = eig(CYY);

      % Multiply left singular vectors with modes matrix to get new modes.
      modes_old = [ui;vi;rhoi];
      modes_old_mod = modes_old*V;

      % Calculate full state variance vector.
      var_full = (modes_old_mod.^2)*diag(D);

      % Extract variance of u, v, rho from this.
      var_u = var_full(1:size(ui,1))*((nu_x/h)^2); var_full(1:size(ui,1))= [];
      var_v = var_full(1:size(vi,1))*((nu_x/h)^2); var_full(1:size(vi,1))= [];
      var_rho = var_full*(DelS*DelS);

      % Plot Everything.

      figure;
      clf;
      subplot(311)
      box on;
      pcolor(XU*h/1e3,YU*h/1e3,var_u(Nodeu(2:app.Ny+1, 2:app.Nx)));
      shading interp;
      colorbar
      title('Variance in U','FontSize',14)

      subplot(312)
      pcolor(XV*h/1e3,YV*h/1e3,var_v(Nodev(2:app.Ny, 2:app.Nx+1)));
      shading interp;
      box on;
      colorbar
      title('Variance in V','FontSize',14)

      subplot(313)
      pcolor(XP*h/1e3,YP*h/1e3,var_rho(NodeRho(2:app.Ny+1, 2:app.Nx+1)));
      shading interp;
      box on;
      colorbar
      title('Variance in S','FontSize',14)

      print('-dpng',imagename)
   end;
end;
