matlab_class/gcmfaces_misc/check_budg.m

function [budg]=check_budg(budgIn,kk);
%function [budg]=check_budg(budgIn,kk);
%
%inputs:
%budgIn is a structure containting tend, trU, trV, 
%    trWtop, trWbot (loading sequence reported below)
%kk is a vertical level (0 for all levels)

gcmfaces_global;

if 0;
  dirIn='nctiles_budg/budgHo/';
  budgIn.tend=read_nctiles([dirIn 'tend/tend'],'tend',1);
  budgIn.trU=read_nctiles([dirIn 'trU/trU'],'trU',1);
  budgIn.trV=read_nctiles([dirIn 'trV/trV'],'trV',1);
  if dirIn(end-1)=='o';
    budgIn.trW=read_nctiles([dirIn 'trW/trW'],'trW',1);
    %
    tend(isnan(tend))=0;
    trU(isnan(trU))=0; trV(isnan(trV))=0;
    trWtop(isnan(trWtop))=0; trWbot(isnan(trWbot))=0;
    %
    budgIn.trWtop=budgIn.trW;
    budgIn.trWbot=budgIn.trW(:,:,2:50);
    budgIn.trWbot(:,:,50)=NaN;
  else;
    budgIn.trWtop=read_nctiles([dirIn 'trWtop/trWtop'],'trWtop',1);
    budgIn.trWbot=read_nctiles([dirIn 'trWbot/trWbot'],'trWbot',1);
  end;
end;%if 0;

tend=budgIn.tend;
trU=budgIn.trU;
trV=budgIn.trV;
trWtop=budgIn.trWtop;
trWbot=budgIn.trWbot;

dt=1;
units='1';
%dt=budgIn.t1-budgIn.t0;
%units = budgIn.specs.units;

%define northern hemisphere as domain of integration
nameMask='Northern Hemisphere';
mask=mygrid.mskC.*mk3D(mygrid.YC>0,mygrid.mskC);
if length(size(tend{1}))<3; mask=mask(:,:,1); end;

%focus on one level
if kk>0;
  mask=mask(:,:,kk);
  tend=tend(:,:,kk);
  trU=trU(:,:,kk);
  trV=trV(:,:,kk);
  trWtop=trWtop(:,:,kk);
  trWbot=trWbot(:,:,kk);
end;

%edit plot title accordingly
descr=nameMask;

%compute northern hemisphere integrals
budg.tend=NaN*dt;
budg.hconv=NaN*dt;
budg.zconv=NaN*dt;

%compute flux convergence
hconv=calc_UV_conv(trU,trV,{});
zconv=-(trWtop-trWbot);%note the confusing upward sign convention ...

%compute sum over domain
budg.tend=nansum(tend.*mask);
budg.hconv=nansum(hconv.*mask);
budg.zconv=nansum(zconv.*mask);

budg
budg.tend-budg.hconv-budg.zconv

1	gforget	1.2	function [budg]=check_budg(budgIn,kk);
2			%function [budg]=check_budg(budgIn,kk);
3			%
4			%inputs:
5			%budgIn is a structure containting tend, trU, trV,
6			% trWtop, trWbot (loading sequence reported below)
7			%kk is a vertical level (0 for all levels)
8	gforget	1.1
9			gcmfaces_global;
10
11	gforget	1.2	if 0;
12			dirIn='nctiles_budg/budgHo/';
13			budgIn.tend=read_nctiles([dirIn 'tend/tend'],'tend',1);
14			budgIn.trU=read_nctiles([dirIn 'trU/trU'],'trU',1);
15			budgIn.trV=read_nctiles([dirIn 'trV/trV'],'trV',1);
16			if dirIn(end-1)=='o';
17			budgIn.trW=read_nctiles([dirIn 'trW/trW'],'trW',1);
18			%
19			tend(isnan(tend))=0;
20			trU(isnan(trU))=0; trV(isnan(trV))=0;
21			trWtop(isnan(trWtop))=0; trWbot(isnan(trWbot))=0;
22			%
23			budgIn.trWtop=budgIn.trW;
24			budgIn.trWbot=budgIn.trW(:,:,2:50);
25			budgIn.trWbot(:,:,50)=NaN;
26			else;
27			budgIn.trWtop=read_nctiles([dirIn 'trWtop/trWtop'],'trWtop',1);
28			budgIn.trWbot=read_nctiles([dirIn 'trWbot/trWbot'],'trWbot',1);
29			end;
30			end;%if 0;
31	gforget	1.1
32			tend=budgIn.tend;
33			trU=budgIn.trU;
34			trV=budgIn.trV;
35			trWtop=budgIn.trWtop;
36			trWbot=budgIn.trWbot;
37
38			dt=1;
39	gforget	1.2	units='1';
40			%dt=budgIn.t1-budgIn.t0;
41			%units = budgIn.specs.units;
42	gforget	1.1
43			%define northern hemisphere as domain of integration
44			nameMask='Northern Hemisphere';
45			mask=mygrid.mskC.*mk3D(mygrid.YC>0,mygrid.mskC);
46			if length(size(tend{1}))<3; mask=mask(:,:,1); end;
47
48	gforget	1.2	%focus on one level
49			if kk>0;
50			mask=mask(:,:,kk);
51			tend=tend(:,:,kk);
52			trU=trU(:,:,kk);
53			trV=trV(:,:,kk);
54			trWtop=trWtop(:,:,kk);
55			trWbot=trWbot(:,:,kk);
56			end;
57
58	gforget	1.1	%edit plot title accordingly
59			descr=nameMask;
60
61			%compute northern hemisphere integrals
62			budg.tend=NaN*dt;
63			budg.hconv=NaN*dt;
64			budg.zconv=NaN*dt;
65
66			%compute flux convergence
67			hconv=calc_UV_conv(trU,trV,{});
68	gforget	1.2	zconv=-(trWtop-trWbot);%note the confusing upward sign convention ...
69
70	gforget	1.1	%compute sum over domain
71			budg.tend=nansum(tend.*mask);
72			budg.hconv=nansum(hconv.*mask);
73			budg.zconv=nansum(zconv.*mask);
74
75			budg
76			budg.tend-budg.hconv-budg.zconv
77