function cthresh(cdfin,cdfout,tind,thresh)
%function cthresh(cdfin,cdfout,tind,thresh)
%
% CTHRESH computes the fraction of time during which the concentrations
%    values at each location exceed a specified threshold for a 
%    an ecomsi model run.  cdfin is the ecomsi.cdf run, and
%    cdfout is a dummy ecomsi.cdf run that has at least n time
%    steps of varout that you don't mind writing over, where n is
%    the length of the vector thresh.
%  
%  tind is the indices of time steps to average over (1 is first step)
varin='conc';
varout='conc';
cdfid1=mexcdf('open',cdfin,'nowrite');
cdfid2=mexcdf('open',cdfout,'write');

[nam,nx]=mexcdf('diminq',cdfid1,'xpos');
[nam,ny]=mexcdf('diminq',cdfid1,'ypos');
[nam,nz]=mexcdf('diminq',cdfid1,'zpos');
nn=length(thresh);
for j=[1:nz-1],
  for k=1:nn,
    st(k)=zeros(nx,ny);
  end
  for i=tind,
    s=mexcdf('varget',cdfid1,varin,[i-1 j-1 0 0],[1 1 ny nx]);
    for k=1:nn,
      sind=s>thresh(k);
      st(k)=st(k)+sind;
    end
  end
  for k=1:nn,
    st(k)=100*st(k)/length(tind);
    mexcdf('varput',cdfid2,varout,[0 k-1 0 0],[1 1 ny nx],st(k));
  end
end
mexcdf('close',cdfid1);
mexcdf('close',cdfid2);
