goodends.m

一个研究声多普勒计程仪很好的工具箱

%  data found for each beam above
   minens=min(inens);   % First good ensemble to use
   maxens=max(oxens);   % Last good ensemble to use
   disp(['Based on the # of good bins, ' num2str(minens) ' is the first good ensemble']);
   disp(['and ' num2str(maxens) ' is the last good, out of ' num2str(nens) ' total ensembles'])
   disp(['from the file ' theMaskFile]);
   disp([' '])
   close(f)

% Now check if Earth coordinates are used.  
% Earth coordinates will not have a masked file.
case 'EARTH'
   for k = 1:4
      v{k} = g{['vel' int2str(k)]};  % Pulls out the velocity data and puts into cell array
   end
   
   
   CHUNK = 100;     % For use below, to break up processing into 'chunks' of 100 for easier/quicker processing
   inens = ones(4,1);   % For use below, to store the first good ensemble for each beam
   oxens = ones(4,1);   % For use below, to store the last good ensemble for each beam
   
   for k = 1:4
      theVarname = ['vel' int2str(k)];   % Pulls out each beam of velocity data
      disp([' checking ' theVarname '...']);
      p = v{k};  % Creates new ncvariable for each beam of data
      m = size(p,1);  % Finds length of ncvariable (i.e. ensemble length)
      m = m(1);
      i = 0;
      t = zeros(m,1);  % Creates a temp array of zeros for good/bad data flags for each beam and ensemble to be stored
      fV = fillval(p);  %Finds the fill value for velocity data
      
      % Below checks velocity data for bad ensembles by seeing if more than 25% of the bins are fill values.
      % A sum of each beam data for each ensemble (each column) is computed.
      % A 'bad' value is found by finding what the value would be if 25% of the bins were fill values.
      % If the sum of the beam data for an ensemble is greater than the 'bad' value, it is marked 'bad'. 
      % For example, if the column has a length of 92 (meaning 92 depth bins), 25% of the column has fill values, or
      %     the 'bad' value, would be 2.3000e+36.
      % If the sum of the column exceeded 2.300e+36, than the column would be marked bad.
      % This check is done in 100 ensemble increments for easier/quicker processing      
      while i < m
         j = i+1:min(i+CHUNK,nens);  % Check data in 'chunks' of 100 ensembles
         q = t(j,:);
         s = sum(p(j,:),2);  % Sums the column of beam
         bad = s > (fV * nbins) * 0.25;  % if # of bins with fill values is greater than 25, marks with '1'
         q = q | bad;
         t(j,:) = q;  % Fills in the temp array with flags for 'good' and 'bad' data
         i = i + CHUNK;  
      end
      
      id = find(t==0);  % Finds the good velocity data for each beam, for each ensemble
      
      if isempty(id)
         id = nan;
      end
    
      % The first/last good ensemble of velocity data (inens/oxens) for each beam is saved here
      inens(k,:) = id(1);     
      oxens(k,:) = id(end);   
    end  %for k=1:4 loop
  
% Will give conservative estimates for the min and max ensembles using the first/last good ensembles of velocity
%  data found for each beam above
  minens = min(inens);      % First good ensemble to use
  maxens = max(oxens);      % Last good ensemble to use
  disp(['Based on the # of good bins, ' num2str(minens) ' is the first good ensemble']);
  disp(['and ' num2str(maxens) ' is the last good ensemble, out of ' num2str(nens) ' total ensembles'])
  disp(['from the file ' theBeamFile]);
  disp([' '])
end     %for case loop
   
%now let's check based on tilt
pitch=g{'Ptch'}(:);
pitch=pitch(minens:maxens);
roll=g{'Roll'}(:);
roll=roll(minens:maxens);
tnens=length(pitch);

dP=diff(pitch);
dr=diff(roll);

%look for pitch and roll that vary by less than 2 degrees
idP=find(abs(dP) < 2 & abs(dr) < 2);

if isempty(idP)
   minp = 1;
elseif isequal(idP(1),1) & ~isequal(idP(2),2)
   minp = idP(2);
else
   minp = idP(1);
end   


if idP(end) >= tnens
   maxp=tnens;
elseif idP(end) < maxens
   maxp=idP(end);
else
   maxp=NaN;
end

minens=minens+minp-1;
maxens=maxp;
   disp(['Based on pitch and roll, ' num2str(minens) ' is the first good ensemble']);
   disp(['and ' num2str(maxens) ' is the last good, out of ' num2str(nens) ' total ensembles'])
   disp(['in file ' theBeamFile])
   disp([' '])


%Needs recently created Netcdf file to run
%some attributes may not be available in older files   
%let's check these results by getting time in gregorian
   gtD=gregorian(time(minens));
   minT=datestr(datenum(gtD(1),gtD(2),gtD(3)),1)
   gtR=gregorian(time(maxens));
   maxT=datestr(datenum(gtR(1),gtR(2),gtR(3)),1)
      
   if ~isempty([ddate,rdate]);
      % If these are not deployment or recovery dates check the less conservative min and maxes
      if ~isequal(datenum(minT),datenum(ddate))
         for ii=2:4;
            inens=sort(inens);
            NgtD=gregorian(time(inens(ii)));
            NminT=datestr(datenum(NgtD(1),NgtD(2),NgtD(3)),1);
            
            if isequal(datenum(minT),datenum(ddate)) 
               disp(['The first good ensemble' num2str(minens) ' recorded on ' num2str(minT)])
               disp(['based on number of good bins, roll, and pitch matches the deployment date ' ])
               minens=inens(ii);
               minT=NminT;			
               break            
               
            elseif isequal(ii,4)
               disp('The first good ensemble was chosen based on deployment date')
               if ~isempty(ddens)
            		minens = ddens;
         		end
					minT = ddate;
				%else
             %  disp('Could not match deployment date')
            end          
         end %for ii
      end
      
      if ~isequal(datenum(maxT),datenum(rdate))   
         for ii=1:3;
            oxens=sort(oxens);            
            NgtR=gregorian(time(oxens(ii)));
            NmaxT=datestr(datenum(NgtR(1),NgtR(2),NgtR(3)),1);
            
            if isequal(datenum(NmaxT),datenum(rdate)) 
               disp('The last good ensemble was chosen based on recovery date')
               disp(['Number of good bins, roll, and pitch gave ensemble ' num2str(minens) ' recorded on ' minT])
               maxens=oxens(ii);
               maxT=NmaxT;            
               break            
            %else
             %  disp('Could not match recovery date')
            end
         end     
      else
         disp(['The last good ensemble is same as recovery date, ' num2str(maxT)]);
      end
      
      disp(' ')
      disp(['The first good ensemble occurs at ' minT])
      disp(['and the recorded deployment date is ' ddate]);
      
      if ~isequal(lower(maxT),lower(rdate))
         %if running batch default to recovery date
         if running(batch)
            gooddate = rdate;
         else
            End_ensemble_dates.Choose_Date = {[maxT ' For Tilt'],[rdate ' Recovery']};
	         Ens_dates = uigetinfo(End_ensemble_dates,'End_ensemble_dates.mat','Date_of_last_good_ensemble');
   	      goodinfo = getinfo(Ens_dates,'Choose_Date');
      	   gooddate = goodinfo(1:11);
         end
         
      switch lower(gooddate)
      case {lower(maxT)}
         maxens = maxens;	
         disp([' '])   
      	disp(['The last good ensemble occurs at ' maxT])
      	disp(['and the recorded recovery date is ' rdate]);

      case {lower(rdate)}
         if ~isempty(rdens)
            maxens = rdens;
         end
         disp([' '])   
         disp(['The recorded recovery date ' rdate])
         disp('was used to determine the last good ensemble')

      otherwise
         disp('invalid recovery date') 
      end
      end
      
   end
   
   close(g)   
   
   %let's cut down the record based on the minens and maxens given
   if ~isequal(1,minens) | ~isequal(maxens,nens)
   	disp(['Trimming file from ensemble ' num2str(minens) ' to ' num2str(maxens)]);
		nctrim(theBeamFile, trimFile, (minens:maxens), 'ensemble','isVerbose');
      
      thecomment=sprintf...
         ('%s\n','Ensembles recorded pre and post deployment were trimmed by goodends.m.');
   	history(trimFile,thecomment);
	else
      disp('no trimming was required')
      if isunix
      	eval(['!cp  ' theBeamFile ' ' trimFile])
		elseif any(findstr(lower(computer), 'pcwin')) | isVMS
      	eval(['!copy  ' theBeamFile ' ' trimFile])
		elseif any(findstr(lower(computer), 'mac')) & ...
			exist('aduplicate') == 2
			feval('aduplicate', theBeamFile, trimFile)
		else
			fcopy(theBeamFile, trimFile)
      end

      thecomment=sprintf('%s\n','No ensembles were trimmed by goodends.m;');
   	history(trimFile,thecomment);
   end

ncclose