Back to home page

MITgcm

 
 

    


Warning, /verification/natl_box/matlab/comp_output.m is written in an unsupported language. File is not indexed.

view on githubraw file Latest commit 2b811582 on 2000-11-14 03:57:32 UTC
cdf4d1244c Patr*0001 % Compare output of new and reference (c32) codes.
                0002 
                0003 p1='../output/';                     % reference (c32)output location
2b8115827f Patr*0004 p2='../../../exe/';                  % new output file location
cdf4d1244c Patr*0005 
                0006 lat=13:2:43; lon=322:2:360;          % latitude, longitude
                0007 
                0008 % model depths
                0009 dZ=[10 10 15 20 20 25 35 50 75 100 150 200 275 ...
                0010       350 415 450 500 500 500 500 500 500 500];
                0011 dpt=dZ/2;
                0012 for i=2:length(dZ)
                0013   dpt(i)=dpt(i)+sum(dZ(1:(i-1)));
                0014 end
                0015 
                0016 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
                0017 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
                0018 
                0019 % load and compare temperature
                0020 
                0021 tx=4;                                % time index
                0022 dp=1;                                % depth level
                0023 cx=[15 28]; cxd=[-1 1]*.005;         % color axes
                0024 
                0025 T1=readbin([p1 'T.001.001.data'],[20 16 23 4],1);
                0026 T2=readbin([p2 'T.001.001.data'],[20 16 23 4],1);
                0027 
                0028 figure(1), clf reset
                0029 set(gcf,'PaperOrientation','portrait')
                0030 set(gcf,'PaperPosition',[0.5 0.5 7.5 10.])
                0031 subplot(311), pcolor(lon,lat,T1(:,:,dp,tx)');
                0032 shading flat, caxis(cx), colorbar
                0033 title(['c32 temperature (deg C), day 30' ...
                0034       ', ' int2str(dpt(dp)) ' m depth' ...
                0035       ', ' ' min=' num2str(min(min(T1(:,:,dp,tx))),4) ...
                0036       ', ' ' max=' num2str(max(max(T1(:,:,dp,tx))),4) ])
                0037 ylabel('Latitude North')
                0038 
                0039 subplot(312), pcolor(lon,lat,T2(:,:,dp,tx)');
                0040 shading flat, caxis(cx), colorbar
                0041 title(['new temperature (deg C), day 30' ...
                0042       ', ' ' min=' num2str(min(min(T1(:,:,dp,tx))),4) ...
                0043       ', ' ' max=' num2str(max(max(T1(:,:,dp,tx))),4) ])
                0044 ylabel('Latitude North')
                0045 
                0046 subplot(313), pcolor(lon,lat,T2(:,:,dp,tx)'-T1(:,:,dp,tx)');
                0047 shading flat, caxis(cxd), colorbar
                0048 title(['Difference (deg C)' ...
                0049       ', ' ' min=' num2str(min(min(T2(:,:,dp,tx)'-T1(:,:,dp,tx)')),4) ...
                0050       ', ' ' max=' num2str(max(max(T2(:,:,dp,tx)'-T1(:,:,dp,tx)')),4) ])
                0051 ylabel('Latitude North'), xlabel('Longitude East')
                0052 orient tall
                0053 filename = 'comp_temp.eps'
                0054 eval([ 'print -depsc ', filename ])
                0055 
                0056 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
                0057 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
                0058 
                0059 % load and compare boundary layer depth
                0060 
                0061 tx=4;                                % time index
                0062 cx=[0 55]; cxd=[-1 1]*4;             % color axes
                0063 
                0064 H1=readbin([p1 'KPPhbl.001.001.data'],[20 16 1 4],1);
                0065 H2=readbin([p2 'KPPhbl.001.001.data'],[20 16 1 4],1);
                0066 
                0067 figure(2), clf reset
                0068 set(gcf,'PaperOrientation','portrait')
                0069 set(gcf,'PaperPosition',[0.5 0.5 7.5 10.])
                0070 subplot(311), pcolor(lon,lat,H1(:,:,1,tx)');
                0071 shading flat, caxis(cx), colorbar
                0072 title(['c32 boundary layer depth (m), day 30' ...
                0073       ', ' ' min=' num2str(min(min(H1(:,:,1,tx))),4) ...
                0074       ', ' ' max=' num2str(max(max(H1(:,:,1,tx))),4) ])
                0075 ylabel('Latitude North')
                0076 
                0077 subplot(312), pcolor(lon,lat,H2(:,:,1,tx)');
                0078 shading flat, caxis(cx), colorbar
                0079 title([ 'new boundary layer depth (m), day 30' ...
                0080       ', ' ' min=' num2str(min(min(H2(:,:,1,tx))),4) ...
                0081       ', ' ' max=' num2str(max(max(H2(:,:,1,tx))),4)] )
                0082 ylabel('Latitude North')
                0083 
                0084 subplot(313), pcolor(lon,lat,H2(:,:,1,tx)'-H1(:,:,1,tx)');
                0085 shading flat, caxis(cxd), colorbar
                0086 title(['Difference (m)' ...
                0087       ', ' ' min=' num2str(min(min(H2(:,:,1,tx)'-H1(:,:,1,tx)')),4) ...
                0088       ', ' ' max=' num2str(max(max(H2(:,:,1,tx)'-H1(:,:,1,tx)')),4)] )
                0089 ylabel('Latitude North'), xlabel('Longitude East')
                0090 orient tall
                0091 filename = 'comp_bldepth.eps'
                0092 eval([ 'print -depsc ', filename ])
                0093 
                0094 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
                0095 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
                0096 
                0097 % load and compare KPP diffusivity
                0098 
                0099 tx=4;                                % time index
                0100 id=2:23;                             % depth index
                0101 cx=[-5 -1.4]; cxd=[-1 1]*.003;       % color axes
                0102 yticks=['.00001';'.0001 ';'.001  ';'.01   '];
                0103 dticks=[3000 1000 300 100 30];
                0104 
                0105 D1=readbin([p1 'KPPdiffKzT.001.001.data'],[20 16 23 4],1);
                0106 D2=readbin([p2 'KPPdiffKzT.001.001.data'],[20 16 23 4],1);
                0107 
                0108 figure(3), clf reset
                0109 set(gcf,'PaperOrientation','portrait')
                0110 set(gcf,'PaperPosition',[0.5 0.5 7.5 10.])
                0111 
                0112 tmp=squeeze(D1(10,:,id,tx))'; tmp(find(~tmp))=1e-10; tmp1=log10(tmp);
                0113 subplot(311), pcolor(lat,-log10(dpt(id)),tmp1);
                0114 shading interp, caxis(cx), h=colorbar;
                0115 set(h,'YTick',-5:-2,'YTickLabel',yticks)
                0116 set(gca,'YTick',-log10(dticks),'YTickLabel',dticks)
                0117 title(['c32 KPP diffusivity (m^2/s), day 30' ...
                0118       ', ' ' min=' num2str(min(min(D1(10,:,id,tx))),4) ...
                0119       ', ' ' max=' num2str(max(max(D1(10,:,id,tx))),4)])
                0120 ylabel('Depth (m)')
                0121 
                0122 tmp=squeeze(D2(10,:,id,tx))'; tmp(find(~tmp))=1e-10; tmp2=log10(tmp);
                0123 subplot(312), pcolor(lat,-log10(dpt(id)),tmp2);
                0124 shading interp, caxis(cx), h=colorbar;
                0125 set(h,'YTick',-5:-2,'YTickLabel',yticks)
                0126 set(gca,'YTick',-log10(dticks),'YTickLabel',dticks)
                0127 title([ 'new KPP diffusivity (m^2/s), day 30' ...
                0128       ', ' ' min=' num2str(min(min(D2(10,:,id,tx))),4) ...
                0129       ', ' ' max=' num2str(max(max(D2(10,:,id,tx))),4) ])
                0130 ylabel('Depth (m)')
                0131 
                0132 subplot(313), pcolor(lat,-log10(dpt(id)),10.^tmp2-10.^tmp1);
                0133 shading interp, caxis(cxd), colorbar
                0134 set(gca,'YTick',-log10(dticks),'YTickLabel',dticks)
                0135 title([ 'Difference (m^2/s)' ...
                0136       ', ' ' min=' num2str(min(min(10.^tmp2-10.^tmp1)),4) ...
                0137       ', ' ' max=' num2str(max(max(10.^tmp2-10.^tmp1)),4) ])
                0138 ylabel('Depth (m)'), xlabel('Latitude North')
                0139 orient tall
                0140 filename = 'comp_diffus.eps'
                0141 eval([ 'print -depsc ', filename ])