monod_1d_carbon/diags/plot_2d.m

clear all


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%% MAKE SOME CHOICES %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
dtype=1;  %1=1d-k; 1.5=1.5d-j; 2=2d-jk
kn=1;       % for 1.5d,2d chose the level toplot
typlot=[1 2 3 7 8 9]; %1=t/par; 2=nutrients, 3=phy/zoo, 4=dom, 5=pom, 6=Rstar
                      % 7=carbon; 8=o2 and alk; 9=surface flux,ph, 10=Chl
dplot=2;  % 1 or 2d

geider=0;

rnp=16;
rsip=16;
rfep=1e-3;
rcp=120;
rop=170;

%%%%%% chose the run directory/time
rundir=['../run/'];
outputdir=['dar1d_20110407_0002/'];

timestr=['0000000000'];

npmax=9;
nzmax=2;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%% read the model output 
ncload([rundir,outputdir,'dic_tave.',timestr,'.t001.nc']);
ncload([rundir,outputdir,'ptr_tave.',timestr,'.t001.nc']);
ncload([rundir,outputdir,'dar_tave.',timestr,'.t001.nc']);
ncload([rundir,outputdir,'tave.',timestr,'.t001.nc']);
%%%%%%%%%%%%%%%%%%%

if (geider==0),
eval(['load ',rundir,'plankton_ini_char_nohead.dat']);
plankton=plankton_ini_char_nohead;
else
eval(['load ',rundir,'gplankton_ini_char_nohead.dat']);
plankton=gplankton_ini_char_nohead;
end


if (dtype==1);   yl=Z;  cay1=200; cay2=0;     end
if (dtype==1.5); xl=Y;  cax1=1; cax2=65;       end


%%%%%%%%%%%%%%%%% MAIN TIME LOOP %%%%%%%%%%%%%%%%%%%%%%
if (dplot==2),
 for itype=typlot,
%%%%%%%%%%
  ncon=5;
  if (itype==1)
% plot T/PAR
   spmax=2;
   tmp1=Ttave; min1=4; max1=12; tstr1='Theta (C)';  np1=1;
   tmp2=dar_PAR_ave; min2=0; max2=220;
                                 tstr2='PAR (\mu Ein/m2/s)'; np2=1;
   ncon=10;
  end % if itype 1
  if (itype==2)
% plot nutrients
   spmax=4;
   tmp1=PO4; min1=0; max1=2; tstr1='PO4 (\muM)';  np1=1;
   tmp2=NO3; min2=0; max2=20; tstr2='NO3 (\muM)';  np2=1;
   tmp3=FeT*1e3; min3=0; max3=1.4; tstr3='Fe (nM)';  np3=1;
   tmp4=SiO2; min4=0; max4=10; tstr4='Si (\muM)';  np4=1;
    ncon=5;
  end % if itype 2
  if (itype==3)
% plot plankton
   spmax=3;
   min1=0; max1=.1; tstr1='Phyto biomass (\muM P)';  np1=npmax;
   min3=0; max3=1e-2; tstr3='Zoo biomass (\muM P)';  np3=nzmax;
   min2=0; max2=1; tstr2='Chl (mg/m3)'; np2=npmax;
  end % if itype 3
  if (itype==4)
% plot DOM
   spmax=3;
   tmp1=DOP; min1=0; max1=0.1; tstr1='DOP (\muM)';  np1=1;
   tmp2=DON; min2=0; max2=3; tstr2='DON (\muM)';  np2=1;
   tmp3=DOFe*1e3; min3=0; max3=0.4; tstr3='DOFe (nM)';  np3=1;
  end % if itype 4
  if (itype==5)
% plot POM
   spmax=4;
   tmp1=POP; min1=0; max1=0.1; tstr1='POP (\muM)';  np1=1;
   tmp2=PON; min2=0; max2=.2; tstr2='PON (\muM)';  np2=1;
   tmp3=POFe*1e3; min3=0; max3=0.1; tstr3='POFe (nM)';  np3=1;
   tmp4=POSi; min4=0; max4=2; tstr4='POSi (\muM)';  np4=1;
  end % if itype 3
%%%%%%%%%
  if (itype==4)
% plot DOM
   spmax=3;
   tmp1=DOP; min1=0; max1=0.1; tstr1='DOP (\muM)';  np1=1;
   tmp2=DON; min2=0; max2=3; tstr2='DON (\muM)';  np2=1;
   tmp3=DOFe*1e3; min3=0; max3=0.4; tstr3='DOFe (nM)';  np3=1;
  end % if itype 4
  if (itype==5)
% plot plankton
   spmax=4;
   tmp1=POP; min1=0; max1=0.1; tstr1='POP (\muM)';  np1=1;
   tmp2=PON; min2=0; max2=.2; tstr2='PON (\muM)';  np2=1;
   tmp3=POFe*1e3; min3=0; max3=0.1; tstr3='POFe (nM)';  np3=1;
   tmp4=POSi; min4=0; max4=2; tstr4='POSi (\muM)';  np4=1;
  end % if itype 3
%%%%%%%%%
  if (itype==6),
   spmax=2;
   min1=0; max1=.05; tstr1='Rstar (\muM P)';  np1=npmax;
   tmp2=PO4; min2=0; max2=0.05; tstr2='PO4 (\muM)';  np2=1;
  end % if itype 6
%%%
  if (itype==7)
% plot carbon
   spmax=4;
   tmp1=DIC; min1=1800; max1=2200; tstr1='DIC (\muM)';  np1=1;
   tmp2=DOC; min2=0; max2=50; tstr2='DOC (\muM)';  np2=1;
   tmp3=POC; min3=0; max3=5; tstr3='POC (nM)';  np3=1;
   tmp4=PIC; min4=0; max4=2; tstr4='PIC (\muM)';  np4=1;
  end % if itype 3
  if (itype==8)
% plot o2 alk
   spmax=2;
   tmp1=ALK; min1=2250; max1=2350; tstr1='ALK (\muM eq)';  np1=1;
   tmp2=O2; min2=0; max2=500; tstr2='O2 (\muM)';  np2=1;
  end % if itype 6
  if (itype==9)
% plot fluxes
   spmax=4;
   tmp1=dic_SUR_ave; min1=-4; max1=4; tstr1='CO2 Flux (mmol/m2/s)';  np1=1;
   tmp2=dic_SURO_ave; min2=-4; max2=4; tstr2='O2 Flux (mmol/,2/s)';  np2=1;
   tmp3=dic_pH_ave; min1=7.8; max1=8.2; tstr3='pH';  np1=1;
   tmp4=dic_pCO2_ave; min2=260e-6; max2=300e-6; tstr4='pCO2';  np2=1;
  end % if itype 9
% plot Chl
  if (itype==10)
% plot chl
   spmax=1;
   tmp1=dar_Chl_ave; min1=0; max1=1; tstr1='Chl (mg/m3)';  np1=1;
  end % if itype 10


%% Plot time/depth output

    figure(itype)
    clf

    if (itype~=3 & itype~=9) % loop through figures
        for isp=1:spmax % loop through subplots
            subplot(spmax,1,isp)

            eval(['plottmp=tmp' num2str(isp) ';']);
%           contourf(T./31104000,[0;-Z(1:end-1)],plottmp',ncon),axis ij
            pcolor(T./31104000,[0;-Z(1:end-1)],plottmp'),axis ij
            shading flat, 
            ylim([0 cay1])

            eval(['title(tstr' num2str(isp) ')'])
            eval(['caxis([min' num2str(isp) ',max' num2str(isp) '])'])
            colorbar
            colormap(flipud(hot))
        end
    end
   if itype==3 % loop through figures
        for isp=1:spmax % loop through subplots
            eval(['npn=np',num2str(isp),';']);


            if isp==1; ptype='Phy'; tag=''; end
            if isp==3; ptype='ZOO'; tag='P'; end
%           if isp==2; ptype='Chl'; tag=''; end
            for ip=1:npn % sum all phyto/zoo-plankton types
                if ip==1;
                    eval(['planktmp=zeros(size(' ptype '1' tag '));']);
                end
                eval(['planktmp=planktmp+' ptype num2str(ip) tag ';']);
            end
            subplot(spmax,1,isp)

            %contourf(T./31104000,[0;-Z(1:end-1)],planktmp',ncon),axis ij
            pcolor(T./31104000,[0;-Z(1:end-1)],planktmp'),axis ij
            shading flat, 
            ylim([0 cay1])

            eval(['title(tstr' num2str(isp) ')'])
            eval(['caxis([min' num2str(isp) ',max' num2str(isp) '])'])
            colorbar
            colormap(flipud(hot))
        end
    end
    if (itype==9) % loop through figures
        for isp=1:spmax % loop through subplots
            subplot(spmax,1,isp)

            eval(['plottmp=tmp' num2str(isp) ';']);
%           contourf(T./31104000,[0;-Z(1:end-1)],plottmp',ncon),axis ij
            plot(T./31104000,plottmp','k-'),

            eval(['title(tstr' num2str(isp) ')'])
%           eval(['caxis([min' num2str(isp) ',max' num2str(isp) '])'])
        end
    end

end

end % if dplot
1	jahn	1.1	clear all
2
3
4			%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
5			%%%% MAKE SOME CHOICES %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
6			dtype=1; %1=1d-k; 1.5=1.5d-j; 2=2d-jk
7			kn=1; % for 1.5d,2d chose the level toplot
8			typlot=[1 2 3 7 8 9]; %1=t/par; 2=nutrients, 3=phy/zoo, 4=dom, 5=pom, 6=Rstar
9			% 7=carbon; 8=o2 and alk; 9=surface flux,ph, 10=Chl
10			dplot=2; % 1 or 2d
11
12			geider=0;
13
14			rnp=16;
15			rsip=16;
16			rfep=1e-3;
17			rcp=120;
18			rop=170;
19
20			%%%%%% chose the run directory/time
21			rundir=['../run/'];
22			outputdir=['dar1d_20110407_0002/'];
23
24			timestr=['0000000000'];
25
26			npmax=9;
27			nzmax=2;
28			%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
29
30			%% read the model output
31			ncload([rundir,outputdir,'dic_tave.',timestr,'.t001.nc']);
32			ncload([rundir,outputdir,'ptr_tave.',timestr,'.t001.nc']);
33			ncload([rundir,outputdir,'dar_tave.',timestr,'.t001.nc']);
34			ncload([rundir,outputdir,'tave.',timestr,'.t001.nc']);
35			%%%%%%%%%%%%%%%%%%%
36
37			if (geider==0),
38			eval(['load ',rundir,'plankton_ini_char_nohead.dat']);
39			plankton=plankton_ini_char_nohead;
40			else
41			eval(['load ',rundir,'gplankton_ini_char_nohead.dat']);
42			plankton=gplankton_ini_char_nohead;
43			end
44
45
46			if (dtype==1); yl=Z; cay1=200; cay2=0; end
47			if (dtype==1.5); xl=Y; cax1=1; cax2=65; end
48
49
50			%%%%%%%%%%%%%%%%% MAIN TIME LOOP %%%%%%%%%%%%%%%%%%%%%%
51			if (dplot==2),
52			for itype=typlot,
53			%%%%%%%%%%
54			ncon=5;
55			if (itype==1)
56			% plot T/PAR
57			spmax=2;
58			tmp1=Ttave; min1=4; max1=12; tstr1='Theta (C)'; np1=1;
59			tmp2=dar_PAR_ave; min2=0; max2=220;
60			tstr2='PAR (\mu Ein/m2/s)'; np2=1;
61			ncon=10;
62			end % if itype 1
63			if (itype==2)
64			% plot nutrients
65			spmax=4;
66			tmp1=PO4; min1=0; max1=2; tstr1='PO4 (\muM)'; np1=1;
67			tmp2=NO3; min2=0; max2=20; tstr2='NO3 (\muM)'; np2=1;
68			tmp3=FeT*1e3; min3=0; max3=1.4; tstr3='Fe (nM)'; np3=1;
69			tmp4=SiO2; min4=0; max4=10; tstr4='Si (\muM)'; np4=1;
70			ncon=5;
71			end % if itype 2
72			if (itype==3)
73			% plot plankton
74			spmax=3;
75			min1=0; max1=.1; tstr1='Phyto biomass (\muM P)'; np1=npmax;
76			min3=0; max3=1e-2; tstr3='Zoo biomass (\muM P)'; np3=nzmax;
77			min2=0; max2=1; tstr2='Chl (mg/m3)'; np2=npmax;
78			end % if itype 3
79			if (itype==4)
80			% plot DOM
81			spmax=3;
82			tmp1=DOP; min1=0; max1=0.1; tstr1='DOP (\muM)'; np1=1;
83			tmp2=DON; min2=0; max2=3; tstr2='DON (\muM)'; np2=1;
84			tmp3=DOFe*1e3; min3=0; max3=0.4; tstr3='DOFe (nM)'; np3=1;
85			end % if itype 4
86			if (itype==5)
87			% plot POM
88			spmax=4;
89			tmp1=POP; min1=0; max1=0.1; tstr1='POP (\muM)'; np1=1;
90			tmp2=PON; min2=0; max2=.2; tstr2='PON (\muM)'; np2=1;
91			tmp3=POFe*1e3; min3=0; max3=0.1; tstr3='POFe (nM)'; np3=1;
92			tmp4=POSi; min4=0; max4=2; tstr4='POSi (\muM)'; np4=1;
93			end % if itype 3
94			%%%%%%%%%
95			if (itype==4)
96			% plot DOM
97			spmax=3;
98			tmp1=DOP; min1=0; max1=0.1; tstr1='DOP (\muM)'; np1=1;
99			tmp2=DON; min2=0; max2=3; tstr2='DON (\muM)'; np2=1;
100			tmp3=DOFe*1e3; min3=0; max3=0.4; tstr3='DOFe (nM)'; np3=1;
101			end % if itype 4
102			if (itype==5)
103			% plot plankton
104			spmax=4;
105			tmp1=POP; min1=0; max1=0.1; tstr1='POP (\muM)'; np1=1;
106			tmp2=PON; min2=0; max2=.2; tstr2='PON (\muM)'; np2=1;
107			tmp3=POFe*1e3; min3=0; max3=0.1; tstr3='POFe (nM)'; np3=1;
108			tmp4=POSi; min4=0; max4=2; tstr4='POSi (\muM)'; np4=1;
109			end % if itype 3
110			%%%%%%%%%
111			if (itype==6),
112			spmax=2;
113			min1=0; max1=.05; tstr1='Rstar (\muM P)'; np1=npmax;
114			tmp2=PO4; min2=0; max2=0.05; tstr2='PO4 (\muM)'; np2=1;
115			end % if itype 6
116			%%%
117			if (itype==7)
118			% plot carbon
119			spmax=4;
120			tmp1=DIC; min1=1800; max1=2200; tstr1='DIC (\muM)'; np1=1;
121			tmp2=DOC; min2=0; max2=50; tstr2='DOC (\muM)'; np2=1;
122			tmp3=POC; min3=0; max3=5; tstr3='POC (nM)'; np3=1;
123			tmp4=PIC; min4=0; max4=2; tstr4='PIC (\muM)'; np4=1;
124			end % if itype 3
125			if (itype==8)
126			% plot o2 alk
127			spmax=2;
128			tmp1=ALK; min1=2250; max1=2350; tstr1='ALK (\muM eq)'; np1=1;
129			tmp2=O2; min2=0; max2=500; tstr2='O2 (\muM)'; np2=1;
130			end % if itype 6
131			if (itype==9)
132			% plot fluxes
133			spmax=4;
134			tmp1=dic_SUR_ave; min1=-4; max1=4; tstr1='CO2 Flux (mmol/m2/s)'; np1=1;
135			tmp2=dic_SURO_ave; min2=-4; max2=4; tstr2='O2 Flux (mmol/,2/s)'; np2=1;
136			tmp3=dic_pH_ave; min1=7.8; max1=8.2; tstr3='pH'; np1=1;
137			tmp4=dic_pCO2_ave; min2=260e-6; max2=300e-6; tstr4='pCO2'; np2=1;
138			end % if itype 9
139			% plot Chl
140			if (itype==10)
141			% plot chl
142			spmax=1;
143			tmp1=dar_Chl_ave; min1=0; max1=1; tstr1='Chl (mg/m3)'; np1=1;
144			end % if itype 10
145
146
147
148
149			%% Plot time/depth output
150
151			figure(itype)
152			clf
153
154			if (itype~=3 & itype~=9) % loop through figures
155			for isp=1:spmax % loop through subplots
156			subplot(spmax,1,isp)
157
158			eval(['plottmp=tmp' num2str(isp) ';']);
159			% contourf(T./31104000,[0;-Z(1:end-1)],plottmp',ncon),axis ij
160			pcolor(T./31104000,[0;-Z(1:end-1)],plottmp'),axis ij
161			shading flat,
162			ylim([0 cay1])
163
164			eval(['title(tstr' num2str(isp) ')'])
165			eval(['caxis([min' num2str(isp) ',max' num2str(isp) '])'])
166			colorbar
167			colormap(flipud(hot))
168			end
169			end
170			if itype==3 % loop through figures
171			for isp=1:spmax % loop through subplots
172			eval(['npn=np',num2str(isp),';']);
173
174
175			if isp==1; ptype='Phy'; tag=''; end
176			if isp==3; ptype='ZOO'; tag='P'; end
177			% if isp==2; ptype='Chl'; tag=''; end
178			for ip=1:npn % sum all phyto/zoo-plankton types
179			if ip==1;
180			eval(['planktmp=zeros(size(' ptype '1' tag '));']);
181			end
182			eval(['planktmp=planktmp+' ptype num2str(ip) tag ';']);
183			end
184			subplot(spmax,1,isp)
185
186			%contourf(T./31104000,[0;-Z(1:end-1)],planktmp',ncon),axis ij
187			pcolor(T./31104000,[0;-Z(1:end-1)],planktmp'),axis ij
188			shading flat,
189			ylim([0 cay1])
190
191			eval(['title(tstr' num2str(isp) ')'])
192			eval(['caxis([min' num2str(isp) ',max' num2str(isp) '])'])
193			colorbar
194			colormap(flipud(hot))
195			end
196			end
197			if (itype==9) % loop through figures
198			for isp=1:spmax % loop through subplots
199			subplot(spmax,1,isp)
200
201			eval(['plottmp=tmp' num2str(isp) ';']);
202			% contourf(T./31104000,[0;-Z(1:end-1)],plottmp',ncon),axis ij
203			plot(T./31104000,plottmp','k-'),
204
205			eval(['title(tstr' num2str(isp) ')'])
206			% eval(['caxis([min' num2str(isp) ',max' num2str(isp) '])'])
207			end
208			end
209
210			end
211
212			end % if dplot