From 7790148cd12947ecf6e16c2fe72f1c7ac4d0e738 Mon Sep 17 00:00:00 2001
From: Antoine Cyril David Hoffmann <ahoffman@spcpc606.epfl.ch>
Date: Tue, 17 May 2022 15:21:42 +0200
Subject: [PATCH] scripts update

---
 matlab/compile_results.m           | 113 +++++++++++----------
 matlab/compute/compute_fa_1D.m     |   4 +-
 matlab/compute/compute_fa_2D.m     |   4 +-
 matlab/compute/mode_growth_meter.m | 155 +++++++++++++++--------------
 matlab/plot/plot_fa.m              |   8 +-
 matlab/plot/plot_fa_gene.m         |  14 +--
 matlab/plot/plot_metric.m          |  26 +++--
 wk/analysis_3D.m                   |  16 +--
 wk/analysis_header.m               |   3 +-
 wk/gene_analysis_3D.m              |  25 +++--
 10 files changed, 196 insertions(+), 172 deletions(-)

diff --git a/matlab/compile_results.m b/matlab/compile_results.m
index d0892abf..952dec78 100644
--- a/matlab/compile_results.m
+++ b/matlab/compile_results.m
@@ -208,63 +208,68 @@ while(CONTINUE)
     JOBNUM   = JOBNUM + 1;
 end
 
-%% Build grids
-
-Nky = numel(ky); 
-if Nky > 1
-    dky = ky(2); 
-    Ly = 2*pi/dky;   
+if(JOBFOUND == 0)
+    disp('no results found, please verify the paths');
+    return;
 else
-    dky = 0;
-    Ly  = 0;   
-end
-[~,iky0] = min(abs(ky)); 
-Ny = 2*Nky-1;  
-y  = linspace(-Ly/2,Ly/2,Ny+1); y = y(1:end-1);
+    %% Build grids
 
-Nkx = numel(kx);
-if Nkx > 1
-    dkx = kx(2);
-    Lx = 2*pi/dkx;
-else
-    dkx = 0;
-    Lx  = 0;
-end    
-[~,ikx0] = min(abs(kx));
-Nx = Nkx;      
-x  = linspace(-Lx/2,Lx/2,Nx+1); x = x(1:end-1);
+    Nky = numel(ky); 
+    if Nky > 1
+        dky = ky(2); 
+        Ly = 2*pi/dky;   
+    else
+        dky = 0;
+        Ly  = 0;   
+    end
+    [~,iky0] = min(abs(ky)); 
+    Ny = 2*Nky-1;  
+    y  = linspace(-Ly/2,Ly/2,Ny+1); y = y(1:end-1);
 
-Nz = numel(z);
+    Nkx = numel(kx);
+    if Nkx > 1
+        dkx = kx(2);
+        Lx = 2*pi/dkx;
+    else
+        dkx = 0;
+        Lx  = 0;
+    end    
+    [~,ikx0] = min(abs(kx));
+    Nx = Nkx;      
+    x  = linspace(-Lx/2,Lx/2,Nx+1); x = x(1:end-1);
 
-[KX,KY] = meshgrid(kx,ky);
-KPERP2 = KX.^2+KY.^2;
-%% Add everything in output structure
-% scaling
-DATA.scale = (1/Nx/Ny)^2;
-% Fields
-DATA.GGAMMA_RI = GGAMMA_; DATA.PGAMMA_RI = PGAMMA_; DATA.HFLUX_X = HFLUX_;
-DATA.Nipj = Nipj_; DATA.Ni00 = Ni00_; DATA.DENS_I = DENS_I_; DATA.TEMP_I = TEMP_I_;
-if(KIN_E)
-DATA.Nepj = Nepj_; DATA.Ne00 = Ne00_; DATA.DENS_E = DENS_E_; DATA.TEMP_E = TEMP_E_;
-end
-DATA.Ts5D = Ts5D_; DATA.Ts3D = Ts3D_; DATA.Ts0D = Ts0D_;
-DATA.PHI  = PHI_; 
-% grids
-DATA.Pe = Pe; DATA.Pi = Pi; 
-DATA.Je = Je; DATA.Ji = Ji; 
-DATA.kx = kx; DATA.ky = ky; DATA.z = z;
-DATA.x  = x;  DATA.y  = y;
-DATA.ikx0 = ikx0; DATA.iky0 = iky0;
-DATA.Nx = Nx; DATA.Ny = Ny; DATA.Nz = Nz; DATA.Nkx = Nkx; DATA.Nky = Nky; 
-DATA.Pmaxe = numel(Pe); DATA.Pmaxi = numel(Pi); DATA.Jmaxe = numel(Je); DATA.Jmaxi = numel(Ji);
-DATA.dir      = DIRECTORY;
-DATA.localdir = DIRECTORY;
-DATA.param_title=['$\nu_{',DATA.CONAME,'}=$', num2str(DATA.NU), ...
-    ', $\kappa_N=$',num2str(DATA.K_N),', $L=',num2str(DATA.L),'$, $N=',...
-    num2str(DATA.Nx),'$, $(P,J)=(',num2str(DATA.PMAXI),',',...
-    num2str(DATA.JMAXI),')$,',' $\mu_{hd}=$(',num2str(DATA.MUx),...
-    ',',num2str(DATA.MUy),')'];
-JOBNUM = LASTJOB;
+    Nz = numel(z);
 
-filename = sprintf([DIRECTORY,'outputs_%.2d.h5'],JOBNUM);
+    [KX,KY] = meshgrid(kx,ky);
+    KPERP2 = KX.^2+KY.^2;
+    %% Add everything in output structure
+    % scaling
+    DATA.scale = (1/Nx/Ny)^2;
+    % Fields
+    DATA.GGAMMA_RI = GGAMMA_; DATA.PGAMMA_RI = PGAMMA_; DATA.HFLUX_X = HFLUX_;
+    DATA.Nipj = Nipj_; DATA.Ni00 = Ni00_; DATA.DENS_I = DENS_I_; DATA.TEMP_I = TEMP_I_;
+    if(KIN_E)
+    DATA.Nepj = Nepj_; DATA.Ne00 = Ne00_; DATA.DENS_E = DENS_E_; DATA.TEMP_E = TEMP_E_;
+    end
+    DATA.Ts5D = Ts5D_; DATA.Ts3D = Ts3D_; DATA.Ts0D = Ts0D_;
+    DATA.PHI  = PHI_; 
+    % grids
+    DATA.Pe = Pe; DATA.Pi = Pi; 
+    DATA.Je = Je; DATA.Ji = Ji; 
+    DATA.kx = kx; DATA.ky = ky; DATA.z = z;
+    DATA.x  = x;  DATA.y  = y;
+    DATA.ikx0 = ikx0; DATA.iky0 = iky0;
+    DATA.Nx = Nx; DATA.Ny = Ny; DATA.Nz = Nz; DATA.Nkx = Nkx; DATA.Nky = Nky; 
+    DATA.Pmaxe = numel(Pe); DATA.Pmaxi = numel(Pi); DATA.Jmaxe = numel(Je); DATA.Jmaxi = numel(Ji);
+    DATA.dir      = DIRECTORY;
+    DATA.localdir = DIRECTORY;
+    DATA.param_title=['$\nu_{',DATA.CONAME,'}=$', num2str(DATA.NU), ...
+        ', $\kappa_N=$',num2str(DATA.K_N),', $L=',num2str(DATA.L),'$, $N=',...
+        num2str(DATA.Nx),'$, $(P,J)=(',num2str(DATA.PMAXI),',',...
+        num2str(DATA.JMAXI),')$,',' $\mu_{hd}=$(',num2str(DATA.MUx),...
+        ',',num2str(DATA.MUy),')'];
+    JOBNUM = LASTJOB;
+
+    filename = sprintf([DIRECTORY,'outputs_%.2d.h5'],JOBNUM);
+end
 end
\ No newline at end of file
diff --git a/matlab/compute/compute_fa_1D.m b/matlab/compute/compute_fa_1D.m
index 79fe35d8..e4774f19 100644
--- a/matlab/compute/compute_fa_1D.m
+++ b/matlab/compute/compute_fa_1D.m
@@ -29,12 +29,12 @@ switch options.SPECIE
 end
 Np = numel(parray); Nj = numel(jarray);
 
-switch options.Z
+switch options.iz
     case 'avg'
         Napj_     = mean(Napj_,5);
         phi_      = mean(data.PHI,3);
     otherwise
-        [~,iz]    = min(abs(options.Z-data.z)); 
+        iz        = options.iz; 
         Napj_     = Napj_(:,:,:,:,iz,:);
         phi_      = data.PHI(:,:,iz);
 end
diff --git a/matlab/compute/compute_fa_2D.m b/matlab/compute/compute_fa_2D.m
index c66ebba5..29a8b093 100644
--- a/matlab/compute/compute_fa_2D.m
+++ b/matlab/compute/compute_fa_2D.m
@@ -26,12 +26,12 @@ switch options.SPECIE
         jarray    = double(data.Ji);
 end
 
-switch options.Z
+switch options.iz
     case 'avg'
         Napj_     = mean(Napj_,5);
         phi_      = mean(data.PHI,3);
     otherwise
-        iz        = options.Z; 
+        iz        = options.iz; 
         Napj_     = Napj_(:,:,:,:,iz,:);
         phi_      = data.PHI(:,:,iz);
 end
diff --git a/matlab/compute/mode_growth_meter.m b/matlab/compute/mode_growth_meter.m
index d0c1a9ce..b8c229d8 100644
--- a/matlab/compute/mode_growth_meter.m
+++ b/matlab/compute/mode_growth_meter.m
@@ -2,98 +2,99 @@ function [FIGURE] = mode_growth_meter(DATA,OPTIONS)
 
 NORMALIZED = OPTIONS.NORMALIZED;
 Nma   = OPTIONS.NMA; %Number moving average
-t  = OPTIONS.TIME;
+
 iz = OPTIONS.iz;
-[~,ikzf] = max(squeeze(mean(abs(squeeze(DATA.PHI(1,:,1,:))),2)));
+[~,ikzf] = max(squeeze(mean(abs(squeeze(DATA.PHI(1,:,iz,:))),2)));
 
 FRAMES = zeros(size(OPTIONS.TIME));
-
 for i = 1:numel(OPTIONS.TIME)
     [~,FRAMES(i)] =min(abs(OPTIONS.TIME(i)-DATA.Ts3D));
 end
+FRAMES = unique(FRAMES);
+t  = DATA.Ts3D(FRAMES);
 
 FIGURE.fig = figure; set(gcf, 'Position',  [100 100 1200 700])
 FIGURE.FIGNAME = 'mode_growth_meter';
 for i = 1:2
-MODES_SELECTOR = i; %(1:Zonal, 2: NZonal, 3: ky=kx)
-
-if MODES_SELECTOR == 1
-    if NORMALIZED
-        plt = @(x,ik) movmean(abs(squeeze(x(1,ik,iz,FRAMES)))./max(abs(squeeze(x(1,ik,iz,FRAMES)))),Nma);
-    else
-        plt = @(x,ik) movmean(abs(squeeze(x(1,ik,iz,FRAMES))),Nma);
-    end
-    kstr = 'k_x';
-    k = DATA.kx;
-    MODESTR = 'Zonal modes';
-elseif MODES_SELECTOR == 2
-    if NORMALIZED
-        plt = @(x,ik) movmean(abs(squeeze(x(ik,1,iz,FRAMES)))./max(abs(squeeze(x(ik,1,iz,FRAMES)))),Nma);
-    else
-        plt = @(x,ik) movmean(abs(squeeze(x(ik,1,iz,FRAMES))),Nma);
-    end
-    kstr = 'k_y';
-    k = DATA.ky;
-    MODESTR = 'NZ modes';
-elseif MODES_SELECTOR == 3
-    if NORMALIZED
-        plt = @(x,ik) movmean(abs(squeeze(x(ik,ik,iz,FRAMES)))./max(abs(squeeze(x(ik,ik,iz,FRAMES)))),Nma);
-    else
-        plt = @(x,ik) movmean(abs(squeeze(x(ik,ik,iz,FRAMES))),Nma);
-    end
-    kstr = 'k_y=k_x';
-    k = DATA.ky;
-    MODESTR = 'Diag modes';
-end 
+    MODES_SELECTOR = i; %(2:Zonal, 1: NZonal, 3: ky=kx)
 
-MODES = 1:numel(k);
-% MODES = zeros(size(OPTIONS.K2PLOT));
-% for i = 1:numel(OPTIONS.K2PLOT)
-%     [~,MODES(i)] =min(abs(OPTIONS.K2PLOT(i)-k));
-% end
+    if MODES_SELECTOR == 2
+        if NORMALIZED
+            plt = @(x,ik) movmean(abs(squeeze(x(1,ik,iz,FRAMES)))./max(abs(squeeze(x(1,ik,iz,FRAMES)))),Nma);
+        else
+            plt = @(x,ik) movmean(abs(squeeze(x(1,ik,iz,FRAMES))),Nma);
+        end
+        kstr = 'k_x';
+        k = DATA.kx;
+        MODESTR = 'Zonal modes';
+    elseif MODES_SELECTOR == 1
+        if NORMALIZED
+            plt = @(x,ik) movmean(abs(squeeze(x(ik,1,iz,FRAMES)))./max(abs(squeeze(x(ik,1,iz,FRAMES)))),Nma);
+        else
+            plt = @(x,ik) movmean(abs(squeeze(x(ik,1,iz,FRAMES))),Nma);
+        end
+        kstr = 'k_y';
+        k = DATA.ky;
+        MODESTR = 'NZ modes';
+    elseif MODES_SELECTOR == 3
+        if NORMALIZED
+            plt = @(x,ik) movmean(abs(squeeze(x(ik,ik,iz,FRAMES)))./max(abs(squeeze(x(ik,ik,iz,FRAMES)))),Nma);
+        else
+            plt = @(x,ik) movmean(abs(squeeze(x(ik,ik,iz,FRAMES))),Nma);
+        end
+        kstr = 'k_y=k_x';
+        k = DATA.ky;
+        MODESTR = 'Diag modes';
+    end 
 
+    MODES = 1:numel(k);
+    % MODES = zeros(size(OPTIONS.K2PLOT));
+    % for i = 1:numel(OPTIONS.K2PLOT)
+    %     [~,MODES(i)] =min(abs(OPTIONS.K2PLOT(i)-k));
+    % end
 
-% plt = @(x,ik) abs(squeeze(x(1,ik,iz,FRAMES)))./max(abs(squeeze(x(1,ik,iz,FRAMES))));
 
-gamma = MODES;
-amp   = MODES;
-i_=1;
-for ik = MODES
-    gr = polyfit(t',log(plt(DATA.PHI,ik)),1);
-    gamma(i_) = gr(1);
-    amp(i_)   = gr(2);
-    i_=i_+1;
-end
+    % plt = @(x,ik) abs(squeeze(x(1,ik,iz,FRAMES)))./max(abs(squeeze(x(1,ik,iz,FRAMES))));
 
-%plot
-subplot(2,3,1+3*(i-1))
-    [YY,XX] = meshgrid(t,fftshift(k,1));
-    pclr = pcolor(XX,YY,abs(plt(fftshift(DATA.PHI,MODES_SELECTOR),1:numel(k))));set(pclr, 'edgecolor','none');  hold on;
-    set(gca,'YDir','normal')
-%     xlim([t(1) t(end)]); %ylim([1e-5 1])
-    xlabel(['$',kstr,'\rho_s$']); ylabel('$t c_s /\rho_s$');
-    title(MODESTR)  
-    
-subplot(2,3,2+3*(i-1))
-    mod2plot = [2:OPTIONS.NMODES+1];
-    for i_ = mod2plot
-        semilogy(t,plt(DATA.PHI,MODES(i_))); hold on;
-%         semilogy(t,exp(gamma(i_).*t+amp(i_)),'--k')
-    end
-    if MODES_SELECTOR == 1
-        semilogy(t,plt(DATA.PHI,ikzf),'--k');
-    end
-    xlim([t(1) t(end)]); %ylim([1e-5 1])
-    xlabel('$t c_s /\rho_s$'); ylabel(['$|\phi_{',kstr,'}|$']);
-    title('Measure time window')
-    
-subplot(2,3,3+3*(i-1))
-    plot(k(MODES),gamma); hold on;
-    plot(k(MODES(mod2plot)),gamma(mod2plot),'x')
-    if MODES_SELECTOR == 1
-        plot(k(ikzf),gamma(ikzf),'ok');
+    gamma = MODES;
+    amp   = MODES;
+    i_=1;
+    for ik = MODES
+        gr = polyfit(t,log(plt(DATA.PHI,ik)),1);
+        gamma(i_) = gr(1);
+        amp(i_)   = gr(2);
+        i_=i_+1;
     end
-    xlabel(['$',kstr,'\rho_s$']); ylabel('$\gamma$');
-    title('Growth rates')
+
+    %plot
+    subplot(2,3,1+3*(i-1))
+        [YY,XX] = meshgrid(t,fftshift(k));
+        pclr = pcolor(XX,YY,abs(plt(fftshift(DATA.PHI,MODES_SELECTOR),1:numel(k))));set(pclr, 'edgecolor','none');  hold on;
+        set(gca,'YDir','normal')
+    %     xlim([t(1) t(end)]); %ylim([1e-5 1])
+        xlabel(['$',kstr,'\rho_s$']); ylabel('$t c_s /\rho_s$');
+        title(MODESTR)  
+
+    subplot(2,3,2+3*(i-1))
+        mod2plot = [2:OPTIONS.NMODES+1];
+        for i_ = mod2plot
+            semilogy(t,plt(DATA.PHI,MODES(i_))); hold on;
+    %         semilogy(t,exp(gamma(i_).*t+amp(i_)),'--k')
+        end
+        if MODES_SELECTOR == 2
+            semilogy(t,plt(DATA.PHI,ikzf),'--k');
+        end
+        xlim([t(1) t(end)]); %ylim([1e-5 1])
+        xlabel('$t c_s /\rho_s$'); ylabel(['$|\phi_{',kstr,'}|$']);
+        title('Measure time window')
+
+    subplot(2,3,3+3*(i-1))
+        plot(k(MODES),gamma); hold on;
+        plot(k(MODES(mod2plot)),gamma(mod2plot),'x')
+        if MODES_SELECTOR == 2
+            plot(k(ikzf),gamma(ikzf),'ok');
+        end
+        xlabel(['$',kstr,'\rho_s$']); ylabel('$\gamma$');
+        title('Growth rates')
 end
 end
\ No newline at end of file
diff --git a/matlab/plot/plot_fa.m b/matlab/plot/plot_fa.m
index 00a2ca2f..efdbd15f 100644
--- a/matlab/plot/plot_fa.m
+++ b/matlab/plot/plot_fa.m
@@ -1,11 +1,11 @@
 function [ FIGURE ] = plot_fa( DATA, OPTIONS )
 
 FIGURE.fig = figure; FIGURE.FIGNAME = ['f_a_',DATA.PARAMS];
-switch OPTIONS.Z
+switch OPTIONS.iz
     case 'avg'
         zcomp = ' z-avg';
     otherwise
-        zcomp = [' z=',sprintf('%2.2f',DATA.z(OPTIONS.Z))];
+        zcomp = [' z=',sprintf('%2.2f',DATA.z(OPTIONS.iz))];
 end
 if OPTIONS.ONED
     [s,x,fsa,fxa] = compute_fa_1D(DATA, OPTIONS); 
@@ -13,12 +13,12 @@ if OPTIONS.ONED
     subplot(1,2,1)
         plot(s,fsa); hold on
         legend(OPTIONS.SPECIE)
-        xlabel('$v_\parallel, (\mu=0)$'); ylabel('$\langle |f_a|^2\rangle_{xy}^{1/2}$'); 
+        xlabel('$v_\parallel, (\mu=0)$'); ylabel(['$\langle |f_a|^2\rangle_{xy}^{1/2}$, ',zcomp]); 
         title(DATA.param_title); 
     subplot(1,2,2)
         plot(x,fxa); hold on;
         legend(OPTIONS.SPECIE)
-        xlabel('$\mu, (v_\parallel=0)$'); ylabel('$\langle |f_a|^2\rangle_{xy}^{1/2}$'); 
+        xlabel('$\mu, (v_\parallel=0)$'); ylabel(['$\langle |f_a|^2\rangle_{xy}^{1/2}$, ',zcomp]);
         if numel(it) == 1
             title(['t=',num2str(DATA.Ts5D(it))]);
         else
diff --git a/matlab/plot/plot_fa_gene.m b/matlab/plot/plot_fa_gene.m
index 997cfe1f..9792d557 100644
--- a/matlab/plot/plot_fa_gene.m
+++ b/matlab/plot/plot_fa_gene.m
@@ -10,13 +10,13 @@ mu = h5read([folder,file],'/coord/mu'); nmu = numel(mu);
 z  = h5read([folder,file],'/coord/z');
 [XX,SS] = meshgrid(mu,vp);
 
-switch OPTIONS.Z
+switch OPTIONS.iz
     case 'avg'
         zcomp_name = ' z-avg';
         zcomp = @(x) squeeze(mean(x,1));
     otherwise
-        zcomp_name = [' z=',sprintf('%2.2f',z(OPTIONS.Z))];
-        zcomp = @(x) squeeze(x(OPTIONS.Z,:,:));
+        zcomp_name = [' z=',sprintf('%2.2f',z(OPTIONS.iz))];
+        zcomp = @(x) squeeze(x(OPTIONS.iz,:,:));
 end
 
 [~,iv0] = min(abs(vp));
@@ -40,14 +40,14 @@ Gdata = Gdata ./ numel(TIMES);
 if OPTIONS.ONED
     switch specie
         case 'e'
-        FFa    = squeeze(Gdata(1,:,:,2));
+        FFa    = zcomp(Gdata(:,:,:,2));
         FFa    = abs(FFa)./max(max(abs(FFa)));
     subplot(1,2,1)
         plot(vp,FFa(:,im0)); hold on;
     subplot(1,2,2)
         plot(mu,FFa(iv0,:)); hold on;
         case 'i'
-        FFa   = squeeze(Gdata(1,:,:,1));    
+        FFa   = zcomp(Gdata(:,:,:,1));    
         FFa    = abs(FFa)./max(max(abs(FFa)));
     end
 
@@ -64,7 +64,7 @@ else
 switch specie
     case 'e'
     name  = '$f_e(v_\parallel,\mu_p)$';
-    FFa    = zcomp(squeeze(Gdata));
+    FFa    = zcomp(squeeze(Gdata(:,:,:,2)));
     FFa    = abs(FFa)./max(max(abs(FFa)));
     switch PLT_FCT
         case 'contour'
@@ -76,7 +76,7 @@ switch specie
     end
     case 'i'
     name  = '$f_i(v_\parallel,\mu_p)$';
-    FFa    = zcomp(squeeze(Gdata));
+    FFa    = zcomp(squeeze(Gdata(:,:,:,1)));
     FFa    = abs(FFa)./max(max(abs(FFa)));
     switch PLT_FCT
         case 'contour'
diff --git a/matlab/plot/plot_metric.m b/matlab/plot/plot_metric.m
index b4b9e080..7b95c932 100644
--- a/matlab/plot/plot_metric.m
+++ b/matlab/plot/plot_metric.m
@@ -1,25 +1,29 @@
 function [ fig ] = plot_metric( data )
 
-names = {'$g^{xx}$','$g^{xy}$','$g^{xz}$','$g^{yy}$','$g^{yz}$','$g^{zz}$',...
-         '$B_0$','$\partial_x B_0$','$\partial_y B_0$','$\partial_z B_0$',...
-         '$J$','$R$','$\phi$','$Z$','$\partial_R x$','$\partial_Z x$'};
-data.geo_arrays = load([data.localdir,'geometry.dat']);
+names = {'Jacobian','gradxB','gradyB','gradzB','gradz_coeff',...
+         'gxx','gxy','gyy','gyz','gzz','hatB','hatR','hatZ'};
+geo_arrays = zeros(2,data.Nz,numel(names));
+
+for i_ = 1:numel(names)
+    namae = names{i_};
+    geo_arrays(:,:,i_) = h5read(data.outfilenames{end},['/data/metric/',namae])';
+end
 fig = figure;
 subplot(311)
-    for i = 1:6
-    plot(data.z, data.geo_arrays(:,i),'DisplayName',names{i}); hold on;
+    for i = 1:5
+    plot(data.z, geo_arrays(1,:,i),'DisplayName',names{i}); hold on;
     end
     xlim([min(data.z),max(data.z)]); legend('show'); title('MoNoLiT geometry');
-    
+
 subplot(312)
-    for i = 7:10
-    plot(data.z, data.geo_arrays(:,i),'DisplayName',names{i}); hold on;
+    for i = 6:10
+    plot(data.z, geo_arrays(1,:,i),'DisplayName',names{i}); hold on;
     end
     xlim([min(data.z),max(data.z)]); legend('show');
 
 subplot(313)
-    for i = 11:16
-    plot(data.z, data.geo_arrays(:,i),'DisplayName',names{i}); hold on;
+    for i = 11:13
+    plot(data.z, geo_arrays(1,:,i),'DisplayName',names{i}); hold on;
     end
     xlim([min(data.z),max(data.z)]); legend('show');
 end
diff --git a/wk/analysis_3D.m b/wk/analysis_3D.m
index 0b8fef01..e9ab6b17 100644
--- a/wk/analysis_3D.m
+++ b/wk/analysis_3D.m
@@ -10,7 +10,7 @@ system(['mkdir -p ',MISCDIR]);
 CMD = ['rsync ', LOCALDIR,'outputs* ',MISCDIR]; disp(CMD);
 system(CMD);
 % Load outputs from jobnummin up to jobnummax
-JOBNUMMIN = 00; JOBNUMMAX = 10;
+JOBNUMMIN = 00; JOBNUMMAX = 08;
 data = compile_results(MISCDIR,JOBNUMMIN,JOBNUMMAX); %Compile the results from first output found to JOBNUMMAX if existing
 
 
@@ -53,7 +53,7 @@ options.PLAN      = 'xy';
 % options.PLAN      = 'sx';
 options.COMP      = 9;
 % options.TIME      = dat.Ts5D;
-options.TIME      = 900:1:990;
+options.TIME      = 920:1:1250;
 data.EPS          = 0.1;
 data.a = data.EPS * 2000;
 create_film(data,options,'.gif')
@@ -99,9 +99,9 @@ if 0
 % options.XPERP     = linspace( 0,6,64);
 options.SPAR      = gene_data.vp';
 options.XPERP     = gene_data.mu';
-options.Z         = 'avg';
-options.T         = 900;
-options.PLT_FCT   = 'contour';
+options.iz        = 13;
+options.T         = 20;
+options.PLT_FCT   = 'pcolor';
 options.ONED      = 0;
 options.non_adiab = 1;
 options.SPECIE    = 'i';
@@ -164,10 +164,10 @@ if 0
 %% Mode evolution
 options.NORMALIZED = 1;
 options.K2PLOT = 1;
-options.TIME   = 5:1:15;
+options.TIME   = 2:20;
 options.NMA    = 1;
-options.NMODES = 5;
-options.iz     = 1;
+options.NMODES = 15;
+options.iz     = 9;
 fig = mode_growth_meter(data,options);
 save_figure(data,fig)
 end
diff --git a/wk/analysis_header.m b/wk/analysis_header.m
index a97d0de7..3ff76bdc 100644
--- a/wk/analysis_header.m
+++ b/wk/analysis_header.m
@@ -6,7 +6,8 @@ outfile ='';
 outfile ='';
 outfile ='';
 % outfile ='shearless_cyclone/128x128x16xdmax_6_L_120_CBC_1.0';
-outfile ='shearless_cyclone/128x128x16xdmax_L_120_CBC_1.0';
+% outfile ='shearless_cyclone/128x128x16xdmax_L_120_CBC_1.0';
+outfile ='shearless_cyclone/linear_CBC_120';
 % outfile ='quick_run/CLOS_1_64x64_5x3_L_120_kN_2.0_kT_0.5_nu_1e-01_SGGK';
 % outfile ='pedestal/64x64x16x2x1_L_300_LnT_20_nu_0.1';
 % outfile ='quick_run/32x32x16_5x3_L_300_q0_2.5_e_0.18_kN_20_kT_20_nu_1e-01_DGGK';
diff --git a/wk/gene_analysis_3D.m b/wk/gene_analysis_3D.m
index 5f83d870..31652e6e 100644
--- a/wk/gene_analysis_3D.m
+++ b/wk/gene_analysis_3D.m
@@ -1,6 +1,7 @@
 % folder = '/misc/gene_results/shearless_cyclone/miller_output_1.0/';
 % folder = '/misc/gene_results/shearless_cyclone/miller_output_0.8/';
-folder = '/misc/gene_results/shearless_cyclone/s_alpha_output_1.0/';
+% folder = '/misc/gene_results/shearless_cyclone/s_alpha_output_1.0/';
+folder = '/misc/gene_results/shearless_cyclone/linear_s_alpha_CBC_100/';
 % folder = '/misc/gene_results/shearless_cyclone/s_alpha_output_0.5/';
 % folder = '/misc/gene_results/shearless_cyclone/LD_s_alpha_output_1.0/';
 % folder = '/misc/gene_results/shearless_cyclone/LD_s_alpha_output_0.8/';
@@ -30,11 +31,11 @@ options.NAME      = '\phi';
 % options.NAME      = 'T_i';
 % options.NAME      = '\Gamma_x';
 % options.NAME      = 'k^2n_e';
-options.PLAN      = 'xy';
+options.PLAN      = 'kxky';
 % options.NAME      ='f_e';
 % options.PLAN      = 'sx';
 options.COMP      = 'avg';
-options.TIME      = [500];
+options.TIME      = [0 500];
 gene_data.a = data.EPS * 2000;
 fig = photomaton(gene_data,options);
 save_figure(gene_data,fig)
@@ -50,7 +51,7 @@ options.NAME      = '\phi';
 % options.NAME      = 'n_i^{NZ}';
 % options.NAME      = '\Gamma_x';
 % options.NAME      = 'n_i';
-options.PLAN      = 'kxky';
+options.PLAN      = 'xy';
 % options.NAME      = 'f_e';
 % options.PLAN      = 'sx';
 options.COMP      = 'avg';
@@ -91,9 +92,21 @@ options.times   = 200:600;
 options.specie  = 'i';
 options.PLT_FCT = 'pcolor';
 options.folder  = folder;
-options.Z       = 'avg';
+options.iz      = 1;
 options.FIELD   = '<f_>';
-options.ONED    = 1;
+options.ONED    = 0;
 % options.FIELD   = 'Q_es';
 plot_fa_gene(options);
 end
+
+if 0
+%% Mode evolution
+options.NORMALIZED = 1;
+options.K2PLOT = 1;
+options.TIME   = 50:150;
+options.NMA    = 1;
+options.NMODES = 5;
+options.iz     = 9;
+fig = mode_growth_meter(gene_data,options);
+save_figure(data,fig)
+end
\ No newline at end of file
-- 
GitLab