clean matlab scripts

nml/parameter.nml

@@ -37,7 +37,7 @@ disp	=0
 nbin	=100
 dx	=5000
 Cfl	=1
-name_sim ='temp2'
+name_sim ='test'
 root	= 'output/'
 
 /

output

-/share/work/bauj0001/output
\ No newline at end of file

output/test/parameters.txt

+!
+!		   _______________________
+!		  |                       |
+!		  |	WIM PARAMETERS    |
+!		  |_______________________|
+!______________________________________________________________________________
+! WAVES PARAMETERS:
+!
+! Tm		-> Peak period					[s]
+! Hs		-> Significant wave height			[m]
+! disp		-> Allowing wave dispersion
+!		   0: Wave dispersion is not allowed, 
+!		      group speed is the same at all spectrum 
+!		      frequency (cg=max[cg(w)])
+!		   1: Wave dispersion is allowed
+!------------------------------------------------------------------------------
+&waves_parameters
+
+Tm	=8
+Hs	=1
+disp	=0
+
+/
+!______________________________________________________________________________
+! MODEL PARAMETERS:
+!
+! nbin		-> Number of grid bin
+! dx		-> Spatial resolution				[m]
+! Cfl		-> Courant–Friedrichs–Lewy condition (0<Cfl<1)
+!		   Only in the case where disp=0. The CFL condition
+!		   is needed to calculate the time step.
+! name_sim	-> name of the output file
+! root		-> destination folder for the output file
+!------------------------------------------------------------------------------
+&model_parameter
+
+nbin	=100
+dx	=5000
+Cfl	=1
+name_sim ='test'
+root	= 'output/'
+
+/
+!______________________________________________________________________________
+! SPECTRUM PARAMETERS:
+!
+! init_spec	-> method to build the wave spectrum
+!	   	   2: Swell
+!		   1: JONSWAP spectrum
+!		   0: Bretschneider spectrum
+! nfreq		-> number of frequency bin
+! Tmin		-> Minimum period 				[s]
+! Tmax		-> Maximum period				[s]
+! alpha_s	-> parameter for jonswap spectrum (init_spec=1)
+! beta_s	-> parameter for jonswap spectrum (init_spec=1)
+! gamma_s	-> parameter for jonswap spectrum (init_spec=1)
+! swell_T	-> swell period (init_spec=2)
+! swell_Hs	-> swell significant height (init_spec=2)
+!------------------------------------------------------------------------------
+&spectrum_parameters
+
+init_spec =0
+nfreq	=800
+Tmin	=2.5
+Tmax	=20
+alpha_s =0.2044
+beta_s	=1.2500
+gamma_s	=3.3
+swell_T =19
+swell_Hs=0.09
+
+/
+!______________________________________________________________________________
+! ICE PARAMETERS
+!
+! X_ice		-> Distance of the ice edge			[m]
+! c_cice	-> Ice concentration 
+! ice_thick	-> method for compute the ice thickness
+!		   0: constant thickness 
+!		   1: thickness is a function of distance
+!		      from ice edge
+! hice		-> Ice thickness (if ice_thick=0)		[m]
+! hmax		-> Maximum ice thickness (if ice_thick=1)	[m]
+! Xh		-> distance where h=hmax/2 (if ice_thicl=1)	[m]
+!------------------------------------------------------------------------------
+&ice_parameters
+X_ice	=50000
+cice	=1
+ice_thick =0
+hice	=3
+hmax	=3
+Xh	=200000
+strain_crit =3e-5
+P_c 	    =0.55
+
+/
+!________________________________________________________________________________
+! FSD PARAMETERS
+
+! FSD_sheme     -> method for compute <D>
+!                  0: dumont et al (2011)
+!                  1: power law
+!
+! minfloe	-> minimum size floe to build the floe size categories [m]
+! maxfloe	-> maximum size floe to build the floe size categories [m]
+! nbcat		-> number of floe size categories
+!--------------------------------------------------------------------------------
+&fsd_parameters
+FSD_scheme	=1
+minfloe		=5
+maxfloe		=400
+nbcat		=60
+
+/
+!_________________________________________________________________________________
+! IDT PARAMETERS
+
+!IDT_scheme	-> compute the ice thickness distribution	
+			0: no distribution
+			1: distribution (rayleigh)
+
+!mu_IDT		-> parameter for the distribution
+!mincat_h	-> minimum ice thickness category
+!maxcat_h	-> maximum ice thickness category
+!nbcat_h	-> number of ice thickness categories
+!---------------------------------------------------------------------------------
+&idt_parameters
+IDT_scheme	=1
+mu_IDT		=0.5
+mincat_h	=0.1
+maxcat_h	=10
+nbcat_h		=50
+/
+!________________________________________________________________________________
+
+
+
+
+
+
+
+
+

scripts_matlab/FSTD.m

+function FSTD
+evalin('base','clear all')
+evalin('base','close all')
+
+[filename,filepath]= uigetfile('*.*','All files');
+path=[filepath,'/',filename];
+
+x=ncread(path,'x_axis');
+t=ncread(path,'time');
+om=ncread(path,'omega');
+spectre=ncread(path,'Spectrum');
+Dave=ncread(path,'Dave');
+FSD=ncread(path,'Floe size distribution');
+Fsize=ncread(path,'floe size');
+thick=ncread(path,'Ice thickness');
+conc=ncread(path,'Ice concentration');
+Hs=ncread(path,'significant height');
+IDT=ncread(path,'Ice Thickness Distribution');
+hcat=ncread(path,'thickness categories');
+dx=(x(2)-x(1))*1000;
+assignin('base','x',x)
+assignin('base','FSD',FSD)
+assignin('base','IDT',IDT)
+assignin('base','Fsize',Fsize)
+assignin('base','hcat',hcat)
+assignin('base','dx',dx)
+
+fig=figure('Name','Floe size and Thickness distribution','position',[200 200 1200 1000],'Numbertitle','off');
+set(fig,'color',[0.2 0.2 0.2])
+cmap=coolcolor;
+assignin('base','cmap',cmap)
+
+plot1=subplot(2,2,3:4);
+
+[hAx,hLine1,hLine2] =plotyy(x,thick,x,Dave);
+set(hAx,'color','none')
+set(hAx(1),'ycolor','b')
+set(hAx(2),'ycolor','g')
+set(hLine1,'linewidth',2,'color','b')
+set(hLine2,'linewidth',2,'color','g')
+h2=gca;
+assignin('base','h2',h2)
+axesPosition=get(h2,'Position');
+ylabel(hAx(2),'Mean Floe size [m]','interpreter','latex');
+ylabel(hAx(1),'Mean Ice thickness [m]','interpreter','latex');
+xlabel('Distance [km]','interpreter','latex');
+hold on
+xLimit = [min(x) max(x)];
+yWidth = 0.05;
+xOffset = -yWidth*diff(xLimit)/axesPosition(3);
+h1 = axes('Position',axesPosition+yWidth.*[-1 0 1 0],...
+    'Color','none','XColor','k','YColor','r',...
+    'XLim',xLimit+[xOffset 0],...
+    'XTick',[],'XTickLabel',[],'NextPlot','add');
+plot(h1,x,Hs(find(Hs>0)),'r','linewidth',2);
+assignin('base','h1',h1)
+ylabel('Significant height [m]','interpreter','latex');
+hold on
+curs=datacursormode(fig);
+set(curs,'enable','on');
+set(curs,'UpdateFcn', @plotother)
+
+    function txt = plotother(~,event)
+        
+        x=evalin('base','x');
+        dx=evalin('base','dx');
+        h1=evalin('base','h1');
+        IDT=evalin('base','IDT');
+        FSD=evalin('base','FSD');
+        Fsize=evalin('base','Fsize');
+        hcat=evalin('base','hcat');
+        cmap=evalin('base','cmap');
+        
+        pos = get(event,'Position');
+        txt = {['',num2str(pos(1))]};
+        [~,posx]=min(abs(x-pos(1)));
+        
+        if (gca==h1)
+            
+            ise = evalin( 'base', 'exist(''linex'') == 1' );
+            
+            if ise
+                evalin('base', 'set(linex,''visible'',''off'')')
+            end
+            
+            FSD1=reshape(FSD(posx,posx,:,:),length(Fsize),length(hcat));
+            IDT1=IDT(posx,:);
+            for i=1:length(hcat)
+                FSD1(:,i)=FSD1(:,i)*IDT1(i);
+                nfloe(:,i)=((FSD1(:,i)*(dx)^2)./Fsize.^2) ;
+            end
+            nfloe=nfloe/sum(sum(nfloe));
+            
+            
+            subplot(2,2,1)
+            title('FSTD')
+            colormap(cmap)
+            pcolor(hcat,Fsize,nfloe);
+            shading interp
+            c=colorbar;
+            xlabel('Floe thickness [m]','interpreter','latex')
+            ylabel('Floe size','interpreter','latex')
+            ylabel(c,'normalized number of floes','interpreter','latex')
+            
+            
+            subplot(2,2,2)
+            title('FSD')
+            FSD2(posx,:)=sum(FSD1,2);
+            nfloe2=((sum(FSD1,2)*dx^2)./Fsize.^2)/sum((sum(FSD1,2)*dx^2)./Fsize.^2);
+            plot(Fsize,nfloe2,'linewidth',2,'color','m')
+            set(gca,'color','none')
+            xlabel('Floe size [m]','interpreter','latex')
+            ylabel('normalized number of floes','interpreter','latex')
+            grid on
+            
+            get(h1,'ylim')
+            linex=plot(h1,[pos(1) pos(1)],get(h1,'ylim'),'linewidth',1.5,'color',[0.2 0.2 0.6],'linestyle','--');
+            assignin('base','linex',linex)
+            
+        end
+        
+    end
+
+end
+
+
+
+            
+            
+            
+            
+

src/initialization.f90

@@ -167,7 +167,7 @@ h_sign(1,1)=Hs
         param1=.false.
         param2=.false.
         
-           open(20,file='output/kice.dat')
+           open(20,file='input/kice.dat')
         
         do iii=1,nbcat_h
         do  ii=1,nfreq